Add tradeoff argument to expm_multiply and invm_multiply both in matfaust and pyfaust.

wrapper/matlab/+matfaust/+poly/expm_multiply.m

@@ -7,6 +7,7 @@
 %> @param B the matrix  or vector to be multiplied by the matrix exponential of A.
 %> @param t (real array) time points.
 %> @param 'K', integer (default value is 10) the greatest polynomial degree of the Chebyshev polynomial basis. The greater it is, the better is the approximate accuracy but note that a larger K increases the computational cost.
+%> @param 'tradeoff', str (optional): 'memory' or 'time' to specify what matters the most: a small memory footprint or a small time of execution. It changes the implementation of pyfaust.poly.poly used behind. It can help when the memory size is limited relatively to the value of rel_err or the size of A and B.
 %> @param 'dev', str (optional): the computation device ('cpu' or 'gpu').
 %>
 %>
@@ -33,9 +34,9 @@ function C = expm_multiply(A, B, t, varargin)
 	dev = 'cpu';
 	rel_err = 1e-6; %TODO: to use later when capable to compute K depending on rel_err
 	K = 10;
-	poly_meth = 1;
 	argc = length(varargin);
 	dev = 'cpu';
+	tradeoff = 'time';
 	if(argc > 0)
 		for i=1:2:argc
 			if(argc > i)
@@ -61,11 +62,11 @@ function C = expm_multiply(A, B, t, varargin)
 					else
 						K = tmparg;
 					end
-				case 'poly_meth'
-					if(argc == i || ~ isscalar(tmparg) || ~isreal(tmparg) || tmparg < 0 || tmparg-floor(tmparg) > 0)
-						error('poly_meth argument must be followed by an integer')
+				case 'tradeoff'
+					if(argc == i || ~ isstr(tmparg))
+						error('tradeoff argument must be followed by ''memory'' or ''time''')
 					else
-						poly_meth = tmparg;
+						tradeoff = tmparg;
 					end
 				otherwise
 					if((isstr(varargin{i}) || ischar(varargin{i}))  && ~ strcmp(tmparg, 'cpu') && ~ startsWith(tmparg, 'gpu'))
@@ -74,6 +75,11 @@ function C = expm_multiply(A, B, t, varargin)
 			end
 		end
 	end
+	if (strcmp(tradeoff, 'memory'))
+		poly_meth = 1;
+	else
+		poly_meth = 2; % tradeoff is time
+	end
 	if (~ issparse(A))
 		error('A must be a sparse matrix.')
 	end

wrapper/matlab/+matfaust/+poly/invm_multiply.m

@@ -5,11 +5,12 @@
 %>
 %> @param A the operator whose inverse is of interest (must be a symmetric positive definite sparse matrix).
 %> @param B the dense matrix or vector to be multiplied by the matrix inverse of A.
-%> @param 'rel_err', rel_err (optional): the targeted relative error between the approximate of the action and the action itself (if you were to compute it with inv(A)*x).
+%> @param 'rel_err', real (optional): the targeted relative error between the approximate of the action and the action itself (if you were to compute it with inv(A)*x).
+%> @param 'tradeoff', str (optional): 'memory' or 'time' to specify what matters the most: a small memory footprint or a small time of execution. It changes the implementation of pyfaust.poly.poly used behind. It can help when the memory size is limited relatively to the value of rel_err or the size of A and B.
 %> @param 'max_K', int (optional): the maximum degree of Chebyshev polynomial to use (useful to limit memory consumption).
 %> @param 'dev', str (optional): the device to instantiate the returned Faust ('cpu' or 'gpu').
 %>
-%> @retval AinvB the array which is the approximate action of matrix inverse  of A on B.
+%> @retval AinvB the array which is the approximate action of matrix inverse of A on B.
 %>
 %> @b Example
 %> @code
@@ -31,7 +32,7 @@ function AinvB = invm_multiply(A, B, varargin)
 	dev = 'cpu';
 	rel_err = 1e-6;
 	max_K = inf;
-	poly_meth = 1;
+	tradeoff = 'time';
 	argc = length(varargin);
 	if(argc > 0)
 		for i=1:2:argc
@@ -58,11 +59,11 @@ function AinvB = invm_multiply(A, B, varargin)
 					else
 						max_K = tmparg;
 					end
-				case 'poly_meth'
-					if(argc == i || ~ isscalar(tmparg) || ~isreal(tmparg) || tmparg < 0 || tmparg-floor(tmparg) > 0)
-						error('poly_meth argument must be followed by an integer')
+				case 'tradeoff'
+					if(argc == i || ~ isstr(tmparg))
+						error('tradeoff argument must be followed by ''memory'' or ''time''')
 					else
-						poly_meth = tmparg;
+						tradeoff = tmparg;
 					end
 				otherwise
 					if((isstr(varargin{i}) || ischar(varargin{i}))  && ~ strcmp(tmparg, 'cpu') && ~ startsWith(tmparg, 'gpu'))
@@ -71,6 +72,11 @@ function AinvB = invm_multiply(A, B, varargin)
 			end
 		end
 	end
+	if (strcmp(tradeoff, 'memory'))
+		poly_meth = 1; % tradeoff is memory
+	else
+		poly_meth = 2;
+	end
 	if (~ issparse(A))
 		error('A must be a sparse matrix.')
 	end
@@ -97,7 +103,7 @@ function AinvB = invm_multiply(A, B, varargin)
 	end
 	coeffs(1) = coeffs(1)*.5;
 	if(poly_meth == 2)
-		TB = T*B
+		TB = T*B;
 		AinvB = matfaust.poly.poly(coeffs, TB, 'dev', dev);
 	elseif(poly_meth == 1)
 		AinvB = matfaust.poly.poly(coeffs, T, 'X', B, 'dev', dev);

wrapper/python/pyfaust/poly.py

@@ -155,7 +155,7 @@ def poly(coeffs, basis='chebyshev', L=None, X=None, dev='cpu', impl='native',
             B with X at None).
             dev: the device to instantiate the returned Faust ('cpu' or 'gpu').
             impl: 'native' (by default) for the C++ impl., "py" for the Python impl.
-            out (np.ndarray): if not None the function result is put into this
+            out: (np.ndarray) if not None the function result is put into this
             np.ndarray. Note that out.flags['F_CONTINUOUS'] must be True. Note that this can't work if the function returns a
             Faust.
 
@@ -474,7 +474,7 @@ class FaustPoly(Faust):
         return next(self.gen)
 
 
-def expm_multiply(A, B, t, K=10, dev='cpu', **kwargs):
+def expm_multiply(A, B, t, K=10, tradeoff='time', dev='cpu', **kwargs):
     """
     Computes an approximate of the action of the matrix exponential of A on B using series of Chebyshev polynomials.
 
@@ -494,6 +494,11 @@ def expm_multiply(A, B, t, K=10, dev='cpu', **kwargs):
         K: the greatest polynomial degree of the Chebyshev polynomial basis.
         The greater it is, the better is the approximate accuracy but note that
         a larger K increases the computational cost.
+        tradeoff: 'memory' or 'time' to specify what matters the most: a small
+        memory footprint or a small time of execution. It changes the
+        implementation of pyfaust.poly.poly used behind. It can help when
+        the memory size is limited relatively to the value of K or the size
+        of A and B.
 
     Returns:
         expm_A_B the result of \f$e^{t_k A} B\f$.
@@ -521,7 +526,12 @@ def expm_multiply(A, B, t, K=10, dev='cpu', **kwargs):
     # check A is PSD or at least square
     if A.shape[0] != A.shape[1]:
         raise ValueError('A must be symmetric positive definite.')
-    poly_meth = 1
+    if tradeoff == 'time':
+        poly_meth = 2
+    elif tradeoff == 'memory':
+        poly_meth = 1
+    else:
+        raise ValueError("tradeoff must be 'memory' or 'time'")
     if 'poly_meth' in kwargs:
         if kwargs['poly_meth'] not in [1, 2, '1', '2']:
             raise ValueError('poly_meth must be 1 or 2')
@@ -563,7 +573,7 @@ def expm_multiply(A, B, t, K=10, dev='cpu', **kwargs):
     else:
         return Y
 
-def invm_multiply(A, B, rel_err=1e-6, max_K=np.inf, dev='cpu', **kwargs):
+def invm_multiply(A, B, rel_err=1e-6, tradeoff='time', max_K=np.inf, dev='cpu', **kwargs):
     """
     Computes an approximate of the action of the matrix inverse of A on B using Chebyshev polynomials.
 
@@ -574,8 +584,12 @@ def invm_multiply(A, B, rel_err=1e-6, max_K=np.inf, dev='cpu', **kwargs):
         (ndarray).
         rel_err: the targeted relative error between the approximate of the action and the action itself (if you were to compute it with np.inv(A)@x).
         max_K: (int, optional) the maximum degree of Chebyshev polynomial to use (useful to limit memory consumption).
-        dev: (optional) 'cpu' or 'gpu',  selects the device to use (currently
-        only 'cpu' is supported).
+        dev: (optional) 'cpu' or 'gpu',  selects the device to use.
+        tradeoff: 'memory' or 'time' to specify what matters the most: a small
+        memory footprint or a small time of execution. It changes the
+        implementation of pyfaust.poly.poly used behind. It can help when
+        the memory size is limited relatively to the value of rel_err or the size
+        of A and B.
 
     Returns:
         The np.ndarray which is the approximate action of matrix inverse  of A on B.
@@ -599,7 +613,12 @@ def invm_multiply(A, B, rel_err=1e-6, max_K=np.inf, dev='cpu', **kwargs):
         raise TypeError('A must be a csr_matrix')
     if A.shape[0] != A.shape[1]:
         raise ValueError('A must be symmetric positive definite.')
-    poly_meth = 1
+    if tradeoff == 'time':
+        poly_meth = 2
+    elif tradeoff == 'memory':
+        poly_meth = 1
+    else:
+        raise ValueError("tradeoff must be 'memory' or 'time'")
     if 'poly_meth' in kwargs:
         poly_meth = kwargs['poly_meth']
         if poly_meth not in [1, 2]: