matfaust::poly Namespace Reference

The matfaust module for polynomial basis as Faust objects. More...

Classes
class	FaustPoly
	Subclass of Faust specialized for orthogonal polynomial basis. More...

Functions
function	basis (L, K, basis_name, varargin)
	Builds the Faust of the polynomial basis defined on the sparse matrix L. More...
function	expm_multiply (A, B, t, varargin)
	Computes an approximate of the action of the matrix exponential of A on B using series of Chebyshev polynomials. More...
function	next (F)
	Gives the next Faust basis of dimension (n+1) from the Faust F polynomial basis of dimension n. More...
function	poly (coeffs, basis, varargin)
	Computes the linear combination of the polynomials defined by basis. More...

Detailed Description

The matfaust module for polynomial basis as Faust objects.

Note

This module is still in BETA status.

Function Documentation

basis()

function matfaust::poly::basis	(	L	,
		K	,
		basis_name	,
		varargin
	)

Builds the Faust of the polynomial basis defined on the sparse matrix L.

Parameters

L	the sparse square matrix.
K	the degree of the last polynomial, i.e. the K+1 first polynomials are built.
basis_name	'chebyshev', and others yet to come.
'T0',matrix	(optional): a sparse matrix to replace the identity as a 0-degree polynomial of the basis.
'dev',str	(optional): the computing device ('cpu' or 'gpu').
'dtype',str	(optional): to decide in which data type the resulting Faust will be encoded ('float' or 'double' by default).

Return values

F	the Faust of the basis composed of the K+1 orthogonal polynomials.

Example

% in a matlab terminal
>> import matfaust.poly.*
>> L = sprand(50, 50, .2);
>> L = L*L';
>> K = 3;
>> F = basis(L, K, 'chebyshev')
 
F =
 
Faust size 200x50, density 0.6612, nnz_sum 6612, 4 factor(s):
- FACTOR 0 (double) SPARSE, size 200x150, density 0.0751333, nnz 2254
- FACTOR 1 (double) SPARSE, size 150x100, density 0.146933, nnz 2204
- FACTOR 2 (double) SPARSE, size 100x50, density 0.4208, nnz 2104
- FACTOR 3 (double) SPARSE, size 50x50, density 0.02, nnz 50
 identity matrix flag

By default, the 0-degree polynomial is the identity. However it is possible to replace the corresponding matrix by any sparse matrix T0 of your choice (with the only constraint that size(T0,1) == size(L, 1)). In that purpose, do as follows:

>> rng(42)
>> T0 = sprand(50, 2, .3);
>> F = basis(L, 3, 'chebyshev', 'T0', T0)
F =
 
Faust size 200x2, density 16.53, nnz_sum 6612, 4 factor(s):
- FACTOR 0 (double) SPARSE, size 200x150, density 0.0751333, nnz 2254
- FACTOR 1 (double) SPARSE, size 150x100, density 0.146933, nnz 2204
- FACTOR 2 (double) SPARSE, size 100x50, density 0.4208, nnz 2104
- FACTOR 3 (double) SPARSE, size 50x2, density 0.5, nnz 50

expm_multiply()

function matfaust::poly::expm_multiply	(	A	,
		B	,
		t	,
		varargin
	)

Computes an approximate of the action of the matrix exponential of A on B using series of Chebyshev polynomials.

Parameters

A	the operator whose exponential is of interest (must be a symmetric positive definite sparse matrix).
B	the matrix or vector to be multiplied by the matrix exponential of A.
t	(real array) time points.
'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.
'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.
'dev',str	(optional): the computing device ('cpu' or 'gpu').
'dtype',str	(optional): to decide in which data type the resulting array C will be encoded ('float' or 'double' by default).

Return values

C	the approximate of e^{t_k A} B. C is a tridimensional array of size (sizef(A,1), size(B,2), size(t, 1)), each slice C(:,:,i) is the action of the matrix exponentatial of A on B according to the time point t(i).

Example

% in a matlab terminal
>> rng(42)
>> import matfaust.poly.expm_multiply
>> L = sprand(5, 5, .2);
>> L = L*L';
>> x = rand(size(L,2), 1);
>> t = linspace(-.5, -.1, 3);
>> C = expm_multiply(L, x, t)
 
C(:,:,1) =
 
    0.1796
    0.1706
    0.3698
    0.4223
    0.2662
 
C(:,:,2) =
 
    0.1809
    0.2164
    0.4254
    0.4261
    0.2748
 
C(:,:,3) =
 
    0.1825
    0.2721
    0.4893
    0.4300
    0.2852
 
>> norm(C(:,:,1) - expm(t(1) * L) * x)
 
ans =
 
   1.0194e-15
 
>> norm(C(:,:,2) - expm(t(2) * L) * x)
 
ans =
 
   4.9763e-15
 
>> norm(C(:,:,3) - expm(t(3) * L) * x)
 
ans =
 
   6.4160e-15
 
>>

see also matlab builtin expm

next()

function matfaust::poly::next

(

F

)

Gives the next Faust basis of dimension (n+1) from the Faust F polynomial basis of dimension n.

Parameters

F	The polynomial basis Faust (must have been generated with matfaust.poly.basis). @reval G the basis of dimension (n+1) with one additional factor added to those of F.

Note

The identical factors between F and G are just views of each other (they are not duplicated in memory).

Example

% in a matlab terminal
>> import matfaust.poly.*
>> rng(42)
>> L = sprand(50, 50, .2);
>> L = L*L';
>> K = 3;
>> F = basis(L, K, 'chebyshev');
>> G = next(F);
>> F
 
F =
 
Faust size 200x50, density 0.6624, nnz_sum 6624, 4 factor(s):
- FACTOR 0 (double) SPARSE, size 200x150, density 0.0752667, nnz 2258
- FACTOR 1 (double) SPARSE, size 150x100, density 0.1472, nnz 2208
- FACTOR 2 (double) SPARSE, size 100x50, density 0.4216, nnz 2108
- FACTOR 3 (double) SPARSE, size 50x50, density 0.02, nnz 50
 identity matrix flag
 
>> G
 
G =
 
Faust size 250x50, density 0.71456, nnz_sum 8932, 5 factor(s):
- FACTOR 0 (double) SPARSE, size 250x200, density 0.04616, nnz 2308
- FACTOR 1 (double) SPARSE, size 200x150, density 0.0752667, nnz 2258
- FACTOR 2 (double) SPARSE, size 150x100, density 0.1472, nnz 2208
- FACTOR 3 (double) SPARSE, size 100x50, density 0.4216, nnz 2108
- FACTOR 4 (double) SPARSE, size 50x50, density 0.02, nnz 50
 identity matrix flag

poly()

function matfaust::poly::poly	(	coeffs	,
		basis	,
		varargin
	)

Computes the linear combination of the polynomials defined by basis.

Parameters

coeffs	the linear combination coefficients (vector).
basis	either the name of the polynomial basis to build on L or the basis if already built externally (as a Faust or an equivalent full array).
'L',matrix	the sparse matrix on which the polynomial basis is built if basis is not already a Faust or a full array.
'X',matrix	if X is set, the linear combination of basis*X is computed (note that the memory space is optimized compared to the manual way of doing first B = basis*X and then calling poly on B without X set).
'dev',str	(optional): the computating device ('cpu' or 'gpu').
'dtype',str	(optional): to decide in which data type the resulting Faust or array will be encoded ('float' or 'double' by default). If basis is a Faust or an array its dtype/class is prioritary over this parameter which is in fact useful only if basis is the name of the basis (a str/char array).

Return values

LC	The linear combination Faust or full array depending on if basis is itself a Faust or a np.ndarray.

Example

>> import matfaust.poly.*
>> L = sprand(50, 50, .02);
>> L = L*L';
>> coeffs = [.5 1 2 3];
>> G = poly(coeffs, 'chebyshev', 'L', L)
 
G =
 
Faust size 50x50, density 0.35, nnz_sum 875, 5 factor(s):
- FACTOR 0 (double) SPARSE, size 50x200, density 0.02, nnz 200
- FACTOR 1 (double) SPARSE, size 200x150, density 0.00916667, nnz 275
- FACTOR 2 (double) SPARSE, size 150x100, density 0.015, nnz 225
- FACTOR 3 (double) SPARSE, size 100x50, density 0.025, nnz 125
- FACTOR 4 (double) SPARSE, size 50x50, density 0.02, nnz 50
 identity matrix flag
 
>>

Which is equivalent to do as below (in two times):

>> K = 3;
>> F = basis(L, K, 'chebyshev');
>> G = poly(coeffs, F)
 
G =
 
Faust size 50x50, density 0.35, nnz_sum 875, 5 factor(s):
- FACTOR 0 (double) SPARSE, size 50x200, density 0.02, nnz 200
- FACTOR 1 (double) SPARSE, size 200x150, density 0.00916667, nnz 275
- FACTOR 2 (double) SPARSE, size 150x100, density 0.015, nnz 225
- FACTOR 3 (double) SPARSE, size 100x50, density 0.025, nnz 125
- FACTOR 4 (double) SPARSE, size 50x50, density 0.02, nnz 50
 identity matrix flag
 
>>

Above G is a Faust because F is too. Below the full array of the Faust F is passed, so an array is returned into GA.

>> GA = poly(coeffs, full(F));
>> ismatrix(GA)
 
ans =
 
        logical
 
        1
 
>>

But of course they are equal:

>> norm(GA-full(G))/(norm(GA)) < 1e-16
 
ans =
 
        logical
 
        1
 
>>