Mentions légales du service

Skip to content
Snippets Groups Projects
Commit 10cb1982 authored by hhakim's avatar hhakim
Browse files

Refactor pyfaust/matfaust.FaustFactory into a new module:...

Refactor pyfaust/matfaust.FaustFactory into a new module: pyfaust/matfaust.fact with all implications (tests, doc, doc gen).
parent 31636d3c
Branches
Tags
No related merge requests found
Showing
with 1448 additions and 1116 deletions
......@@ -7,10 +7,10 @@ if(BUILD_DOCUMENTATION)
string(CONCAT DOXYGEN_FILE_PATTERNS "*.cpp *.hpp *.h *.cu *.hu")
endif()
if(BUILD_WRAPPER_MATLAB)
string(CONCAT DOXYGEN_FILE_PATTERNS ${DOXYGEN_FILE_PATTERNS} " Faust.m StoppingCriterion.m ConstraintGeneric.m ConstraintMat.m ConstraintReal.m ConstraintInt.m ConstraintName.m ParamsFact.m ParamsHierarchicalFact.m ParamsPalm4MSA.m FaustFactory.m hadamard.m quickstart.m fft.m bsl.m runtimecmp.m runall.m version.m faust_fact.m ParamsHierarchicalFactSquareMat.m ParamsHierarchicalFactRectMat.m license.m omp.m wht.m dft.m eye.m rand.m")
string(CONCAT DOXYGEN_FILE_PATTERNS ${DOXYGEN_FILE_PATTERNS} " Faust.m StoppingCriterion.m ConstraintGeneric.m ConstraintMat.m ConstraintReal.m ConstraintInt.m ConstraintName.m ParamsFact.m ParamsHierarchicalFact.m ParamsPalm4MSA.m FaustFactory.m hadamard.m quickstart.m fft.m bsl.m runtimecmp.m runall.m version.m faust_fact.m ParamsHierarchicalFactSquareMat.m ParamsHierarchicalFactRectMat.m license.m omp.m wht.m dft.m eye.m rand.m eigtj.m fact_hierarchical.m fact.m fact_palm4msa.m fgft_givens.m fgft_palm.m svdtj.m ")
endif()
if(BUILD_WRAPPER_PYTHON)
string(CONCAT DOXYGEN_FILE_PATTERNS ${DOXYGEN_FILE_PATTERNS} " __init__.py factparams.py demo.py tools.py")
string(CONCAT DOXYGEN_FILE_PATTERNS ${DOXYGEN_FILE_PATTERNS} " __init__.py factparams.py demo.py tools.py fact.py")
endif()
configure_file(${FAUST_DOC_SRC_DIR}/Doxyfile.in ${PROJECT_BINARY_DIR}/doc/Doxyfile @ONLY)
# ./gen_doc/images/* files is duplicated in doc/html/ to call images documentation in the source code with relative path of image's files, from build directory.
......
......@@ -519,7 +519,7 @@ EXCLUDE_SYMLINKS = NO
EXCLUDE_PATTERNS = *faust_MatDense_gpu.h *faust_MatSparse_gpu.h
EXCLUDE_SYMBOLS = mtimes_trans subsasgn get_factor_nonopt reshape ABC ParamsPalm4MSAFGFT UpdateCholesky UpdateCholeskyFull UpdateCholeskySparse greed_omp_chol
EXCLUDE_SYMBOLS = mtimes_trans subsasgn get_factor_nonopt reshape ABC ParamsPalm4MSAFGFT UpdateCholesky UpdateCholeskyFull UpdateCholeskySparse greed_omp_chol fact_hierarchical_constends fact_palm4msa_constends
# The EXAMPLE_PATH tag can be used to specify one or more files or
# directories that contain example code fragments that are included (see
......
#!/bin/bash
@PYTHON3_EXE@ -m doxypypy.doxypypy -a -c $* | @PYTHON3_EXE@ py_filterout_namespace.py pyfaust.__init__. pyfaust.factparams. pyfaust.demo. pyfaust.tools.
@PYTHON3_EXE@ -m doxypypy.doxypypy -a -c $* | @PYTHON3_EXE@ \
py_filterout_namespace.py pyfaust.__init__. pyfaust.factparams. pyfaust.demo. \
pyfaust.tools. pyfaust.fact.
......@@ -19,7 +19,7 @@ classdef FaustFactoryTest < matlab.unittest.TestCase
methods (Test)
function test_fact_palm4msa(this)
disp('Test FaustFactory.fact_palm4msa()')
disp('Test matfaust.fact.fact_palm4msa()')
import matfaust.*
import matfaust.factparams.*
num_facts = 2;
......@@ -38,7 +38,7 @@ classdef FaustFactoryTest < matlab.unittest.TestCase
stop_crit = StoppingCriterion(200);
params = ParamsPalm4MSA(cons, stop_crit, 'is_update_way_R2L', is_update_way_R2L, 'init_lambda', init_lambda, 'step_size', ParamsFact.DEFAULT_STEP_SIZE,...
'constant_step_size', ParamsFact.DEFAULT_CONSTANT_STEP_SIZE);
F = FaustFactory.fact_palm4msa(M, params)
F = matfaust.fact.fact_palm4msa(M, params)
this.verifyEqual(size(F), size(M))
%disp('norm F: ')
%norm(F, 'fro')
......@@ -51,7 +51,7 @@ classdef FaustFactoryTest < matlab.unittest.TestCase
end
function test_fact_palm4msaCplx(this)
disp('Test FaustFactory.fact_palm4msaCplx()')
disp('Test matfaust.fact.fact_palm4msaCplx()')
import matfaust.*
import matfaust.factparams.*
num_facts = 2;
......@@ -69,7 +69,7 @@ classdef FaustFactoryTest < matlab.unittest.TestCase
cons{2} = ConstraintReal(ConstraintName(ConstraintName.NORMCOL), 32, 32, 1.0);
stop_crit = StoppingCriterion(200);
params = ParamsPalm4MSA(cons, stop_crit, 'is_update_way_R2L', is_update_way_R2L, 'init_lambda', init_lambda);
F = FaustFactory.fact_palm4msa(M, params)
F = matfaust.fact.fact_palm4msa(M, params)
this.verifyEqual(size(F), size(M))
%disp('norm F: ')
%norm(F, 'fro')
......@@ -80,7 +80,7 @@ classdef FaustFactoryTest < matlab.unittest.TestCase
end
function test_fact_hierarchical(this)
disp('Test FaustFactory.fact_hierarchical()')
disp('Test matfaust.fact.fact_hierarchical()')
import matfaust.*
import matfaust.factparams.*
num_facts = 4;
......@@ -106,7 +106,7 @@ classdef FaustFactoryTest < matlab.unittest.TestCase
stop_crit = StoppingCriterion(200);
stop_crit2 = StoppingCriterion(200);
params = ParamsHierarchicalFact(fact_cons, res_cons, stop_crit, stop_crit2);
F = FaustFactory.fact_hierarchical(M, params)
F = matfaust.fact.fact_hierarchical(M, params)
this.verifyEqual(size(F), size(M))
%disp('norm F: ')
%norm(F, 'fro')
......@@ -119,7 +119,7 @@ classdef FaustFactoryTest < matlab.unittest.TestCase
end
function test_fact_hierarchicalCplx(this)
disp('Test FaustFactory.fact_hierarchicalCplx()')
disp('Test matfaust.fact.fact_hierarchicalCplx()')
import matfaust.*
import matfaust.factparams.*
num_facts = 4;
......@@ -146,7 +146,7 @@ classdef FaustFactoryTest < matlab.unittest.TestCase
stop_crit2 = StoppingCriterion(200);
params = ParamsHierarchicalFact(fact_cons, res_cons, stop_crit, stop_crit2,...
'init_lambda', init_lambda, 'is_update_way_R2L', is_update_way_R2L);
F = FaustFactory.fact_hierarchical(M, params)
F = matfaust.fact.fact_hierarchical(M, params)
this.verifyEqual(size(F), size(M))
%disp('norm F: ')
%norm(F, 'fro')
......@@ -159,11 +159,11 @@ classdef FaustFactoryTest < matlab.unittest.TestCase
end
function test_fgft_givens(this)
disp('Test FaustFactory.fgft_givens()')
disp('Test matfaust.fact.fgft_givens()')
import matfaust.*
load([this.faust_paths{1} '../../../misc/data/mat/test_GivensDiag_Lap_U_J.mat'])
% Lap and J available
[F,D] = FaustFactory.fgft_givens(Lap, J);%, 0, 'verbosity', 1);
[F,D] = matfaust.fact.fgft_givens(Lap, J);%, 0, 'verbosity', 1);
this.verifyEqual(size(F), size(Lap))
%disp('norm F: ')
%norm(F, 'fro')
......@@ -173,18 +173,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] = FaustFactory.eigtj(Lap, J);%, 0, 'verbosity', 2);
[F2,D2] = matfaust.fact.eigtj(Lap, J);%, 0, 'verbosity', 2);
this.verifyEqual(full(F2),full(F))
this.verifyEqual(D,D2)
end
function test_fgft_givens_parallel(this)
disp('Test FaustFactory.fgft_givens() -- 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] = FaustFactory.fgft_givens(Lap, J, t); %, 'verbosity', 2);
[F,D] = matfaust.fact.fgft_givens(Lap, J, t); %, 'verbosity', 2);
this.verifyEqual(size(F), size(Lap))
%disp('norm F: ')
%norm(F, 'fro')
......@@ -194,13 +194,13 @@ classdef FaustFactoryTest < matlab.unittest.TestCase
% misc/test/src/C++/GivensFGFTParallel.cpp.in
this.verifyEqual(err,0.0410448, 'AbsTol', 0.00001)
% verify it works the same using the eigtj() alias function
[F2,D2] = FaustFactory.eigtj(Lap, J, t); %, 'verbosity', 2);
[F2,D2] = matfaust.fact.eigtj(Lap, J, t); %, 'verbosity', 2);
this.verifyEqual(full(F2),full(F))
this.verifyEqual(D,D2)
end
function test_fgft_palm(this)
disp('Test FaustFactory.fgft_palm()')
disp('Test matfaust.fact.fgft_palm()')
import matfaust.*
import matfaust.factparams.*
num_facts = 4;
......@@ -225,7 +225,7 @@ classdef FaustFactoryTest < matlab.unittest.TestCase
params.init_lambda = 128;
params = ParamsHierarchicalFact(fact_cons, res_cons, stop_crit, stop_crit2, 'is_fact_side_left', params.fact_side == 1, 'is_update_way_R2L', params.update_way == 1, 'init_lambda', params.init_lambda, 'step_size', params.stepsize, 'constant_step_size', false, 'is_verbose', params.verbose ~= 1);
diag_init_D = diag(init_D)
[F,D,lambda] = FaustFactory.fgft_palm(U, Lap, params, diag_init_D)
[F,D,lambda] = matfaust.fact.fgft_palm(U, Lap, params, diag_init_D)
this.verifyEqual(size(F), size(U))
%disp('norm F: ')
%norm(F, 'fro')
......
......@@ -730,8 +730,8 @@ class TestFaustPyCplx(TestFaustPy):
class TestFaustFactory(unittest.TestCase):
def testFactPalm4MSA(self):
print("Test FaustFactory.fact_palm4msa()")
from pyfaust import FaustFactory;
from pyfaust.fact import fact_palm4msa
print("Test pyfaust.fact.fact_palm4msa()")
from pyfaust.factparams import ConstraintReal,\
ConstraintInt, ConstraintName
from pyfaust.factparams import ParamsPalm4MSA, StoppingCriterion
......@@ -743,7 +743,7 @@ class TestFaustFactory(unittest.TestCase):
# init_facts.append(np.eye(32))
#M = np.random.rand(500, 32)
M = \
loadmat(sys.path[-1]+"/../../../misc/data/mat/config_compared_palm2.mat")['data']
loadmat(sys.path[0]+"/../../../misc/data/mat/config_compared_palm2.mat")['data']
# default step_size
cons1 = ConstraintInt(ConstraintName(ConstraintName.SPLIN), 500, 32, 5)
cons2 = ConstraintReal(ConstraintName(ConstraintName.NORMCOL), 32,
......@@ -752,7 +752,7 @@ class TestFaustFactory(unittest.TestCase):
param = ParamsPalm4MSA([cons1, cons2], stop_crit, init_facts=None,
is_update_way_R2L=False, init_lambda=1.0,
is_verbose=False, constant_step_size=False)
F, _lambda= FaustFactory.fact_palm4msa(M, param, ret_lambda=True)
F, _lambda = fact_palm4msa(M, param, ret_lambda=True)
#F.display()
#print("normF", F.norm("fro"))
self.assertEqual(F.shape, M.shape)
......@@ -765,8 +765,8 @@ class TestFaustFactory(unittest.TestCase):
self.assertAlmostEqual(norm(E,"fro")/norm(M,"fro"), 0.270954109668, places=4)
def testFactHierarch(self):
print("Test FaustFactory.fact_hierarchical()")
from pyfaust import FaustFactory
print("Test pyfaust.fact.fact_hierarchical()")
from pyfaust.fact import fact_hierarchical
from pyfaust.factparams import ParamsHierarchicalFact, StoppingCriterion
from pyfaust.factparams import ConstraintReal, ConstraintInt,\
ConstraintName
......@@ -775,7 +775,7 @@ class TestFaustFactory(unittest.TestCase):
init_lambda = 1.0
#M = np.random.rand(500, 32)
M = \
loadmat(sys.path[-1]+"/../../../misc/data/mat/matrix_hierarchical_fact.mat")['matrix']
loadmat(sys.path[0]+"/../../../misc/data/mat/matrix_hierarchical_fact.mat")['matrix']
# default step_size
fact0_cons = ConstraintInt(ConstraintName(ConstraintName.SPLIN), 500, 32, 5)
fact1_cons = ConstraintInt(ConstraintName(ConstraintName.SP), 32, 32, 96)
......@@ -789,7 +789,7 @@ class TestFaustFactory(unittest.TestCase):
[res0_cons, res1_cons, res2_cons],
stop_crit1, stop_crit2,
is_verbose=False)
F = FaustFactory.fact_hierarchical(M, param)
F = fact_hierarchical(M, param)
self.assertEqual(F.shape, M.shape)
E = F.todense()-M
#print("err.:",norm(F.todense(), "fro"), norm(E,"fro"), norm (M,"fro"),
......@@ -799,8 +799,8 @@ class TestFaustFactory(unittest.TestCase):
self.assertAlmostEqual(norm(E,"fro")/norm(M,"fro"), 0.268513, places=5)
def testFactHierarchCplx(self):
print("Test FaustFactory.fact_hierarchicalCplx()")
from pyfaust import FaustFactory
print("Test pyfaust.fact.fact_hierarchicalCplx()")
from pyfaust.fact import fact_hierarchical
from pyfaust.factparams import ParamsHierarchicalFact, StoppingCriterion
from pyfaust.factparams import ConstraintReal,\
ConstraintInt, ConstraintName
......@@ -809,7 +809,7 @@ class TestFaustFactory(unittest.TestCase):
init_lambda = 1.0
#M = np.random.rand(500, 32)
M = \
loadmat(sys.path[-1]+"/../../../misc/data/mat/matrix_hierarchical_fact.mat")['matrix']
loadmat(sys.path[0]+"/../../../misc/data/mat/matrix_hierarchical_fact.mat")['matrix']
M = M + np.complex(0,1)*M
# default step_size
fact0_cons = ConstraintInt(ConstraintName(ConstraintName.SPLIN), 500, 32, 5)
......@@ -824,7 +824,7 @@ class TestFaustFactory(unittest.TestCase):
[res0_cons, res1_cons, res2_cons],
stop_crit1, stop_crit2,
is_verbose=False)
F = FaustFactory.fact_hierarchical(M, param)
F = fact_hierarchical(M, param)
self.assertEqual(F.shape, M.shape)
#print(F.todense())
E = F.todense()-M
......@@ -835,8 +835,8 @@ class TestFaustFactory(unittest.TestCase):
self.assertAlmostEqual(norm(E,"fro")/norm(M,"fro"), 1.0063, places=4)
def testFactPalm4MSACplx(self):
print("Test FaustFactory.fact_palm4msaCplx()")
from pyfaust import FaustFactory
print("Test pyfaust.fact.fact_palm4msaCplx()")
from pyfaust.fact import fact_palm4msa
from pyfaust.factparams import ParamsPalm4MSA,StoppingCriterion
from pyfaust.factparams import ConstraintReal,\
ConstraintInt, ConstraintName
......@@ -848,7 +848,7 @@ class TestFaustFactory(unittest.TestCase):
# init_facts.append(np.eye(32))
#M = np.random.rand(500, 32)
M = \
loadmat(sys.path[-1]+"/../../../misc/data/mat/config_compared_palm2.mat")['data']
loadmat(sys.path[0]+"/../../../misc/data/mat/config_compared_palm2.mat")['data']
M = M + np.complex(0,1)*M
# default step_size
cons1 = ConstraintInt(ConstraintName(ConstraintName.SPLIN), 500, 32, 5)
......@@ -858,7 +858,7 @@ class TestFaustFactory(unittest.TestCase):
param = ParamsPalm4MSA([cons1, cons2], stop_crit, init_facts=None,
is_verbose=False, constant_step_size=False,
is_update_way_R2L=False, init_lambda=1.0)
F,_lambda = FaustFactory.fact_palm4msa(M, param, ret_lambda=True)
F,_lambda = fact_palm4msa(M, param, ret_lambda=True)
#F.display()
#print("normF", F.norm("fro"))
self.assertEqual(F.shape, M.shape)
......@@ -899,13 +899,14 @@ class TestFaustFactory(unittest.TestCase):
dft(n, False).normalize().toarray()))
def testFGFTGivens(self):
print("Test FaustFactory.fgft_givens()")
from pyfaust import FaustFactory as FF
L = loadmat(sys.path[-1]+"/../../../misc/data/mat/test_GivensDiag_Lap_U_J.mat")['Lap']
print("Test fact.fgft_givens()")
print(sys.path)
import pyfaust.fact
L = loadmat(sys.path[0]+"/../../../misc/data/mat/test_GivensDiag_Lap_U_J.mat")['Lap']
L = L.astype(np.float64)
J = \
int(loadmat(sys.path[-1]+"/../../../misc/data/mat/test_GivensDiag_Lap_U_J.mat")['J'])
F, D = FF.fgft_givens(L, J, 0, verbosity=0)
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)
print("Lap norm:", norm(L, 'fro'))
err = norm((F*D.todense())*F.T.todense()-L,"fro")/norm(L,"fro")
print("err: ", err)
......@@ -914,21 +915,21 @@ class TestFaustFactory(unittest.TestCase):
self.assertAlmostEqual(err, 0.0326529, places=7)
def testFGFTGivensParallel(self):
print("Test FaustFactory.fgft_givens() -- parallel")
from pyfaust import FaustFactory as FF
L = loadmat(sys.path[-1]+"/../../../misc/data/mat/test_GivensDiag_Lap_U_J.mat")['Lap']
print("Test pyfaust.fact.fgft_givens() -- parallel")
from pyfaust.fact import eigtj, fgft_givens
L = loadmat(sys.path[0]+"/../../../misc/data/mat/test_GivensDiag_Lap_U_J.mat")['Lap']
L = L.astype(np.float64)
J = \
int(loadmat(sys.path[-1]+"/../../../misc/data/mat/test_GivensDiag_Lap_U_J.mat")['J'])
int(loadmat(sys.path[0]+"/../../../misc/data/mat/test_GivensDiag_Lap_U_J.mat")['J'])
t = int(L.shape[0]/2)
F, D = FF.fgft_givens(L, J, t, verbosity=0)
F, D = fgft_givens(L, J, 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 = FF.eigtj(L, J, t, verbosity=0)
F2, D2 = eigtj(L, J, 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)
......@@ -937,24 +938,23 @@ class TestFaustFactory(unittest.TestCase):
self.assertEqual(err, err2)
def testFactPalm4MSA_fgft(self):
print("Test FaustFactory._fact_palm4msa_fgft()")
print("Test pyfaust.fact._fact_palm4msa_fgft()")
from pyfaust.factparams import ConstraintReal,\
ConstraintInt, ConstraintName, StoppingCriterion
#from pyfaust.factparams import ParamsPalm4MSAFGFT, StoppingCriterion
from numpy import diag,copy
from pyfaust import FaustFactory
from pyfaust.fact import _fact_palm4msa_fgft
from pyfaust.factparams import ParamsPalm4MSAFGFT
from pyfaust import FaustFactory as FF
L = \
loadmat(sys.path[-1]+"/../../../misc/data/mat/ref_test_PALM4SMA_FFT2")['data']
loadmat(sys.path[0]+"/../../../misc/data/mat/ref_test_PALM4SMA_FFT2")['data']
init_D = \
loadmat(sys.path[-1]+"/../../../misc/data/mat/ref_test_PALM4SMA_FFT2")['p_init_D']
loadmat(sys.path[0]+"/../../../misc/data/mat/ref_test_PALM4SMA_FFT2")['p_init_D']
init_D = copy(diag(init_D))
init_facts1 = \
loadmat(sys.path[-1]+"/../../../misc/data/mat/ref_test_PALM4SMA_FFT2")['p_init_facts1']
loadmat(sys.path[0]+"/../../../misc/data/mat/ref_test_PALM4SMA_FFT2")['p_init_facts1']
init_facts2 = \
loadmat(sys.path[-1]+"/../../../misc/data/mat/ref_test_PALM4SMA_FFT2")['p_init_facts2']
loadmat(sys.path[0]+"/../../../misc/data/mat/ref_test_PALM4SMA_FFT2")['p_init_facts2']
init_facts = [ init_facts1, init_facts2 ]
L = L.astype(np.float64)
# other params should be set from file also but do it manually
......@@ -968,7 +968,7 @@ class TestFaustFactory(unittest.TestCase):
init_D=init_D,
is_update_way_R2L=False, init_lambda=128,
is_verbose=True, step_size=1e-6)
F, D, _lambda = FaustFactory._fact_palm4msa_fgft(L, param, ret_lambda=True)
F, D, _lambda = _fact_palm4msa_fgft(L, param, ret_lambda=True)
print("Lap norm:", norm(L, 'fro'))
print("out lambda:", _lambda)
D = diag(D)
......@@ -979,8 +979,8 @@ class TestFaustFactory(unittest.TestCase):
self.assertAlmostEqual(err, 1.39352e-5, places=5)
def testFactHierarchFGFT(self):
print("Test FaustFactory.fgft_palm()")
from pyfaust import FaustFactory
print("Test pyfaust.fact.fgft_palm()")
import pyfaust.fact
from pyfaust.factparams import ParamsHierarchicalFact, StoppingCriterion
from pyfaust.factparams import ConstraintReal, ConstraintInt,\
ConstraintName
......@@ -990,13 +990,13 @@ class TestFaustFactory(unittest.TestCase):
init_lambda = 1.0
#M = np.random.rand(500, 32)
U = \
loadmat(sys.path[-1]+"/../../../misc/data/mat/HierarchicalFactFFT_test_U_L_params.mat")['U']
loadmat(sys.path[0]+"/../../../misc/data/mat/HierarchicalFactFFT_test_U_L_params.mat")['U']
Lap = \
loadmat(sys.path[-1]+"/../../../misc/data/mat/HierarchicalFactFFT_test_U_L_params.mat")['Lap'].astype(np.float)
loadmat(sys.path[0]+"/../../../misc/data/mat/HierarchicalFactFFT_test_U_L_params.mat")['Lap'].astype(np.float)
init_D = \
loadmat(sys.path[-1]+"/../../../misc/data/mat/HierarchicalFactFFT_test_U_L_params.mat")['init_D']
loadmat(sys.path[0]+"/../../../misc/data/mat/HierarchicalFactFFT_test_U_L_params.mat")['init_D']
params_struct = \
loadmat(sys.path[-1]+'/../../../misc/data/mat/HierarchicalFactFFT_test_U_L_params.mat')['params']
loadmat(sys.path[0]+'/../../../misc/data/mat/HierarchicalFactFFT_test_U_L_params.mat')['params']
nfacts = params_struct['nfacts'][0,0][0,0] #useless
niter1 = params_struct['niter1'][0,0][0,0]
niter2 = params_struct['niter2'][0,0][0,0]
......@@ -1030,7 +1030,7 @@ class TestFaustFactory(unittest.TestCase):
constant_step_size=False)
diag_init_D = copy(diag(init_D))
print("norm(init_D):", norm(init_D))
F,D, _lambda = FaustFactory.fgft_palm(U, Lap, param,
F,D, _lambda = pyfaust.fact.fgft_palm(U, Lap, param,
diag_init_D,
ret_lambda=True)
print("out_lambda:", _lambda)
......@@ -1046,7 +1046,8 @@ if __name__ == "__main__":
if(len(sys.argv)> 1):
# argv[1] is for adding a directory in PYTHONPATH
# (to find pyfaust module)
sys.path.append(sys.argv[1])
#sys.path.append(sys.argv[1])
sys.path = [sys.argv[1]]+sys.path # gives priority to pyfaust path
del sys.argv[1] # deleted to avoid interfering with unittest
from pyfaust import Faust
if(len(sys.argv) > 1):
......
function check_fact_mat(funcname, M)
if(~ ismatrix(M) || isscalar(M))
error([funcname,'() 1st argument (M) must be a matrix.'])
end
if(~ isnumeric(M))
error([funcname, '() 1st argument (M) must be real or complex.'])
end
% if(~ isreal(M))
% error([funcname, '() doesn''t yet support complex matrix factorization.'])
% end
end
%==========================================================================================
%> @brief Computes the eigenvalues and the eigenvectors transform (as a Faust object) using the truncated Jacobi algorithm.
%>
%> The eigenvalues and the eigenvectors are approximate. The trade-off between accuracy and sparsity can be set through the parameters J and t.
%>
%>
%> @b Usage
%>
%> &nbsp;&nbsp;&nbsp; @b eigtj(M, J) calls the non-parallel Givens algorithm.<br/>
%> &nbsp;&nbsp;&nbsp; @b eigtj(M, J, 0) or eigtj(M, J, 1) do the same as in previous line.<br/>
%> &nbsp;&nbsp;&nbsp; @b eigtj(M, J, t) calls the parallel Givens algorithm (if t > 1, otherwise it calls basic Givens algorithm)<br/>
%> &nbsp;&nbsp;&nbsp; @b eigtj(M, J, t, 'verbosity', 2) same as above with a level of verbosity of 2 in output. <br/>
%>
%> @param M the matrix to diagonalize. Must be real and symmetric.
%> @param J 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).
%> @param t 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 > size(M,1)/2 won't lead to the desired effect because the maximum number of rotation matrices per factor is anyway size(M,1)/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.
%> @param verbosity the level of verbosity, the greater the value the more info. is displayed.
%>
%> @retval [V,D]
%> - 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 eigenvalues set in D.
%> - D the approximate sparse diagonal matrix of the eigenvalues (in ascendant order along the rows/columns).
%>
%> @b Example
%> @code
%> import matfaust.fact.eigtj
%>
%> % get a Laplacian to diagonalize
%> load('Laplacian_256_community.mat')
%> % do it
%> [Uhat, Dhat] = eigtj(Lap, size(Lap,1)*100, size(Lap, 1)/2, 'verbosity', 2)
%> % Uhat is the Fourier matrix/eigenvectors approximattion as a Faust (200 factors + permutation mat.)
%> % Dhat the eigenvalues diagonal matrix approx.
%> @endcode
%>
%>
%>
%> @b References:
%> - [1] Le Magoarou L., Gribonval R. and Tremblay N., "Approximate fast
%> graph Fourier transforms via multi-layer sparse approximations",
%> IEEE Transactions on Signal and Information Processing
%> over Networks 2018, 4(2), pp 407-420
%>
%> <p> @b See @b also fact.fgft_givens, fact.fgft_palm
%>
%==========================================================================================
function [V,D] = eigtj(M, J, varargin)
[V, D] = matfaust.fact.fgft_givens(M, J, varargin{:});
V = factors(V,1:numfactors(V))
% copy seems unecessary but it's to workaround a bug (temporarily)
end
%> @package matfaust.fact @brief The matfaust factorization module.
%>
%>
%> This module gives access to the main factorization algorithms of
%> FAuST. Those algorithms can factorize a dense matrix to a sparse product
%> (i.e. a Faust object).
%>
%> There are several factorization algorithms.
%>
%> - The first one is PALM4MSA :
%> which stands for Proximal Alternating Linearized Minimization for
%> Multi-layer Sparse Approximation. Note that PALM4MSA is not
%> intended to be used directly. You should rather rely on the second algorithm.
%>
%> - The second one is the Hierarchical Factorization algorithm:
%> this is the central algorithm to factorize a dense matrix to a Faust.
%> It makes iterative use of Palm4MSA to proceed with the factorization of a given
%> dense matrix.
%>
%> - The third group of algorithms is for FGFT computing: fact.fgft_palm fact.fgft_givens fact.eigtj
%>
%>
% ======================================================================
%==========================================================================================
%> @brief Factorizes the matrix M with Hierarchical Factorization using the parameters set in p.
%>
%>
%> @param M the dense matrix to factorize.
%> @param p is a set of factorization parameters. It might be a fully defined instance of parameters (matfaust.factparams.ParamsHierarchicalFact) or a simplified expression which designates a pre-defined parametrization:
%> - 'squaremat' to use pre-defined parameters typically used to factorize a Hadamard square matrix of order a power of two (see matfaust.demo.hadamard).
%> - {'rectmat', j, k, s} to use pre-defined parameters used for instance in factorization of the MEG matrix which is a rectangular matrix of size m*n such that m < n (see matfaust.demo.bsl); j is the number of factors, k the sparsity of the main factor's columns, and s the sparsity of rows for all other factors except the residuum (that is the first factor here because the factorization is made toward the left -- is_side_fact_left == true, cf. matfaust.factparams.ParamsHierarchicalFact).
%> </br>The residuum has a sparsity of P*rho^(num_facts-1). <br/> By default, rho == .8 and P = 1.4. It's possible to set custom values with for example p == { 'rectmat', j, k, s, 'rho', .4, 'P', .7}. <br/>The sparsity is here the number of non-zero elements.
%> @note - The fully defined parameters (ParamsHierarchicalFact instance) used/generated by the function are available in the return result (so one can consult what precisely mean the simplified parameterizations and possibly adjust the attributes to factorize again).
%> @note - This function has its shorthand matfaust.faust_fact(). For convenience you might use it like this:
%> @code
%> import matfaust.*
%> F = faust_fact(M, p) % equiv. to matfaust.fact.fact_hierarchical(M, p)
%> @endcode
%>
%> @retval F The Faust object result of the factorization.
%> @retval [F, lambda, p_obj] = fact_hierarchical(M, p) to optionally get lambda (scale) and the p_obj ParamsHierarchicalFact instance used to factorize.
%>
%> @b Example 1: Fully Defined Parameters for a Random Matrix Factorization
%> @code
%> import matfaust.*
%> import matfaust.factparams.*
%> import matfaust.fact.fact_hierarchical
%> M = rand(500, 32);
%> fact_cons = ConstraintList('splin', 5, 500, 32, 'sp', 96, 32, 32, 'sp', 96, 32, 32);
%> res_cons = ConstraintList('normcol', 1, 32, 32, 'sp', 666, 32, 32, 'sp', 333, 32, 32);
%> stop_crit = StoppingCriterion(200);
%> stop_crit2 = StoppingCriterion(200);
%> params = ParamsHierarchicalFact(fact_cons, res_cons, stop_crit, stop_crit2, 'is_update_way_R2L', false, 'init_lambda', 1.0);
%> F = fact_hierarchical(M, params)
%> @endcode
%> Faust::HierarchicalFact<FPP,DEVICE>::compute_facts : factorization 1/3<br/>
%> Faust::HierarchicalFact<FPP,DEVICE>::compute_facts : factorization 2/3<br/>
%> Faust::HierarchicalFact<FPP,DEVICE>::compute_facts : factorization 3/3<br/>
%>
%> F =
%>
%> Faust size 500x32, density 0.189063, nnz_sum 3025, 4 factor(s):
%> - FACTOR 0 (real) SPARSE, size 500x32, density 0.15625, nnz 2500
%> - FACTOR 1 (real) SPARSE, size 32x32, density 0.09375, nnz 96
%> - FACTOR 2 (real) SPARSE, size 32x32, density 0.09375, nnz 96
%> - FACTOR 3 (real) SPARSE, size 32x32, density 0.325195, nnz 333
%>
%> @b Example 2: Simplified Parameters for Hadamard Factorization
%>@code
%> import matfaust.*
%> import matfaust.fact.fact_hierarchical
%> % generate a Hadamard Faust of size 32x32
%> FH = wh(32);
%> H = full(FH); % the full matrix version
%> % factorize it
%> FH2 = fact_hierarchical(H, 'squaremat');
%> % test the relative error
%> norm(FH-FH2, 'fro')/norm(FH, 'fro') % the result is 1.1015e-16, the factorization is accurate
%>@endcode
%>FH =
%>
%>Faust size 32x32, density 0.3125, nnz_sum 320, 5 factor(s):
%>- FACTOR 0 (real) SPARSE, size 32x32, density 0.0625, nnz 64
%>- FACTOR 1 (real) SPARSE, size 32x32, density 0.0625, nnz 64
%>- FACTOR 2 (real) SPARSE, size 32x32, density 0.0625, nnz 64
%>- FACTOR 3 (real) SPARSE, size 32x32, density 0.0625, nnz 64
%>- FACTOR 4 (real) SPARSE, size 32x32, density 0.0625, nnz 64
%>
%>FH2 =
%>
%>Faust size 32x32, density 0.3125, nnz_sum 320, 5 factor(s):
%>- FACTOR 0 (real) SPARSE, size 32x32, density 0.0625, nnz 64
%>- FACTOR 1 (real) SPARSE, size 32x32, density 0.0625, nnz 64
%>- FACTOR 2 (real) SPARSE, size 32x32, density 0.0625, nnz 64
%>- FACTOR 3 (real) SPARSE, size 32x32, density 0.0625, nnz 64
%>- FACTOR 4 (real) SPARSE, size 32x32, density 0.0625, nnz 64
%>
%>
%>@b Example 3: Simplified Parameters for a Rectangular Matrix Factorization (the BSL demo MEG matrix)
%>
%> @code
%> >> % in a matlab terminal
%> >> import matfaust.*
%> >> import matfaust.fact.fact_hierarchical
%> >> load('matrix_MEG.mat')
%> >> MEG = matrix;
%> >> num_facts = 9;
%> >> k = 10;
%> >> s = 8;
%> >> MEG16 = fact_hierarchical(MEG, {'rectmat', num_facts, k, s})
%> @endcode
%> MEG16 =
%>
%> Faust size 204x8193, density 0.0631655, nnz_sum 105573, 9 factor(s):
%> - FACTOR 0 (real) SPARSE, size 204x204, density 0.293613, nnz 12219
%> - FACTOR 1 (real) SPARSE, size 204x204, density 0.0392157, nnz 1632
%> - FACTOR 2 (real) SPARSE, size 204x204, density 0.0392157, nnz 1632
%> - FACTOR 3 (real) SPARSE, size 204x204, density 0.0392157, nnz 1632
%> - FACTOR 4 (real) SPARSE, size 204x204, density 0.0392157, nnz 1632
%> - FACTOR 5 (real) SPARSE, size 204x204, density 0.0392157, nnz 1632
%> - FACTOR 6 (real) SPARSE, size 204x204, density 0.0392157, nnz 1632
%> - FACTOR 7 (real) SPARSE, size 204x204, density 0.0392157, nnz 1632
%> - FACTOR 8 (real) SPARSE, size 204x8193, density 0.0490196, nnz 81930
%>
%> @code
%> >> % verify the constraint k == 10, on column 5
%> >> fac9 = factors(MEG16,9);
%> >> numel(nonzeros(fac9(:,5)))
%> @endcode
%> ans =
%>
%> 10
%>
%> @code
%> >> % now verify the s constraint is respected on MEG16 factor 2
%> >> numel(nonzeros(factors(MEG16, 2)))/size(MEG16,1)
%> @endcode
%>
%>ans =
%>
%> 8
%> <p> @b See @b also matfaust.faust_fact, factparams.ParamsHierarchicalFact, factparams.ParamsHierarchicalFactSquareMat, factparams.ParamsHierarchicalFactRectMat
%==========================================================================================
function varargout = fact_hierarchical(M, p)
import matfaust.Faust
import matfaust.factparams.*
import matfaust.fact.check_fact_mat
check_fact_mat('matfaust.fact.fact_hierarchical', M)
if(~ isa(p, 'ParamsHierarchicalFact') && ParamsFactFactory.is_a_valid_simplification(p))
p = ParamsFactFactory.createParams(M, p);
end
mex_constraints = cell(2, p.num_facts-1);
if(~ isa(p ,'ParamsHierarchicalFact'))
error('p must be a ParamsHierarchicalFact object.')
end
%mex_fact_constraints = cell(1, p.num_facts-1)
for i=1:p.num_facts-1
cur_cell = cell(1, 4);
cur_cell{1} = p.constraints{i}.name.conv2str();
cur_cell{2} = p.constraints{i}.param;
cur_cell{3} = p.constraints{i}.num_rows;
cur_cell{4} = p.constraints{i}.num_cols;
%mex_fact_constraints{i} = cur_cell;
mex_constraints{1,i} = cur_cell;
end
%mex_residuum_constraints = cell(1, p.num_facts-1)
for i=1:p.num_facts-1
cur_cell = cell(1, 4);
cur_cell{1} = p.constraints{i+p.num_facts-1}.name.conv2str();
cur_cell{2} = p.constraints{i+p.num_facts-1}.param;
cur_cell{3} = p.constraints{i+p.num_facts-1}.num_rows;
cur_cell{4} = p.constraints{i+p.num_facts-1}.num_cols;
%mex_residuum_constraints{i} = cur_cell;
mex_constraints{2,i} = cur_cell;
end
if(~ p.is_mat_consistent(M))
error('M''s number of columns must be consistent with the last residuum constraint defined in p. Likewise its number of rows must be consistent with the first factor constraint defined in p.')
end
% the setters for num_rows/cols verifies consistency with constraints
mex_params = struct('nfacts', p.num_facts, 'cons', {mex_constraints}, 'niter1', p.stop_crits{1}.num_its,'niter2', p.stop_crits{2}.num_its, 'sc_is_criterion_error', p.stop_crits{1}.is_criterion_error, 'sc_error_treshold', p.stop_crits{1}.error_treshold, 'sc_max_num_its', p.stop_crits{1}.max_num_its, 'sc_is_criterion_error2', p.stop_crits{2}.is_criterion_error, 'sc_error_treshold2', p.stop_crits{2}.error_treshold, 'sc_max_num_its2', p.stop_crits{2}.max_num_its, 'nrow', p.data_num_rows, 'ncol', p.data_num_cols, 'fact_side', p.is_fact_side_left, 'update_way', p.is_update_way_R2L, 'verbose', p.is_verbose, 'init_lambda', p.init_lambda);
if(isreal(M))
[lambda, core_obj] = mexHierarchical_factReal(M, mex_params);
else
[lambda, core_obj] = mexHierarchical_factCplx(M, mex_params);
end
F = Faust(core_obj, isreal(M));
varargout = {F, lambda, p};
end
%==========================================================================================
%> @brief Approximates M by A S_1 ... S_n B using fact_hierarchical.
%>
%> @Example
%> @code
%> import matfaust.fact.fact_hierarchical
%> import matfaust.factparams.*
%>
%> p = ParamsHierarchicalFact(…
%> ConstraintList('spcol', 2, 10, 20, 'sp', 30, 10, 10), ConstraintList('sp', 4, 10, 20, 'splin', 5, 10, 10),…
%> StoppingCriterion(50), StoppingCriterion(50),…
%> 'is_fact_side_left', true, 'is_verbose', false…
%> );
%> M = rand(10,10);
%> A = matfaust.rand(10,10);
%> B = matfaust.rand(20, 10);
%> [F, lamdba, ~] = fact_hierarchical_constends(M, p, A, B)
%>
%> assert(norm(A - factors(F,1))/norm(A) <= eps(double(1)))
%> assert(norm(B - factors(F,4))/norm(B) <= eps(double(1)))
%> @endcode
%>
%==========================================================================================
function varargout = fact_hierarchical_constends(M, p, A, B)
import matfaust.factparams.*
import matfaust.fact.fact_hierarchical
if(~ ismatrix(A) || ~ ismatrix(B))
error('A and B must be matrices.')
end
consA = ConstraintList('const', A, size(A, 1), size(A, 2));
consB = ConstraintList('const', B, size(B, 1), size(B, 2));
consts = p.constraints;
nconsts = length(p.constraints);
% consts: factor constraints + residuum constraints
fac_cons = {};
res_cons = {};
for i=1:p.num_facts-1
fac_cons = { fac_cons{:}, consts{i} };
end
for i=p.num_facts:length(consts)
res_cons = { res_cons{:}, consts{i} };
end
assert(length(fac_cons) == length(res_cons))
% add two constants factor constraints for A and B to the old constraints
% According to the factorization direction, switch A and B positions
if(p.is_fact_side_left)
new_fact_cons = { consB.clist{:}, fac_cons{:} };
new_res_cons = { res_cons{:}, consA.clist{:} };
else
new_fact_cons = { consA.clist{:}, fac_cons{:} };
new_res_cons = { res_cons{:}, consB.clist{:} };
end
p = ParamsHierarchicalFact(new_fact_cons, new_res_cons,...
p.stop_crits{1}, p.stop_crits{2}, 'is_update_way_R2L', p.is_update_way_R2L, ...
'init_lambda', p.init_lambda, 'step_size', p.step_size, 'constant_step_size', ...
p.constant_step_size, 'is_verbose', p.is_verbose, 'is_fact_side_left', p.is_fact_side_left);
[F, lambda, p] = fact_hierarchical(M, p);
f1 = factors(F, 1);
f1 = f1 / lambda;
nF = cell(1, numfactors(F));
nF{1} = f1;
for i=2:numfactors(F)
nF{i} = factors(F, i);
end
nF{2} = nF{2}*lambda;
F = matfaust.Faust(nF);
varargout = {F, lambda, p};
end
%==========================================================================================
%> @brief Factorizes the matrix M with Palm4MSA algorithm using the parameters set in p.
%>
%>
%> @param M the dense matrix to factorize.
%> @param p the ParamsPalm4MSA instance to define the algorithm parameters.
%>
%> @retval F the Faust object result of the factorization.
%> @retval [F, lambda] = fact_palm4msa(M, p) to optionally get lambda (scale).
%>
%> @b Example
%>
%> @code
%> import matfaust.*
%> import matfaust.factparams.*
%> import matfaust.fact.fact_palm4msa
%> M = rand(500, 32);
%> cons = ConstraintList('splin', 5, 500, 32, 'normcol', 1, 32, 32);
%> stop_crit = StoppingCriterion(200);
%> params = ParamsPalm4MSA(cons, stop_crit, 'is_update_way_R2L', false, 'init_lambda', 1.0);
%> F = fact_palm4msa(M, params)
%> @endcode
%>
%> F =
%>
%> Faust size 500x32, density 0.22025, nnz_sum 3524, 2 factor(s):
%> - FACTOR 0 (real) SPARSE, size 500x32, density 0.15625, nnz 2500
%> - FACTOR 1 (real) SPARSE, size 32x32, density 1, nnz 1024
%>
%>
%==========================================================================================
function [F,lambda] = fact_palm4msa(M, p)
import matfaust.Faust
import matfaust.fact.check_fact_mat
mex_constraints = cell(1, length(p.constraints));
check_fact_mat('matfaust.fact.fact_palm4msa', M)
if(~ isa(p ,'matfaust.factparams.ParamsPalm4MSA'))
error('p must be a ParamsPalm4MSA object.')
end
for i=1:length(p.constraints)
cur_cell = cell(1, 4);
cur_cell{1} = p.constraints{i}.name.conv2str();
cur_cell{2} = p.constraints{i}.param;
cur_cell{3} = p.constraints{i}.num_rows;
cur_cell{4} = p.constraints{i}.num_cols;
mex_constraints{i} = cur_cell;
end
if(~ p.is_mat_consistent(M))
error('M''s number of columns must be consistent with the last residuum constraint defined in p. Likewise its number of rows must be consistent with the first factor constraint defined in p.')
end
% put mex_constraints in a cell array again because mex eats one level of array
mex_params = struct('data', M, 'nfacts', p.num_facts, 'cons', {mex_constraints}, 'init_facts', {p.init_facts}, 'niter', p.stop_crit.num_its, 'sc_is_criterion_error', p.stop_crit.is_criterion_error, 'sc_error_treshold', p.stop_crit.error_treshold, 'sc_max_num_its', p.stop_crit.max_num_its, 'update_way', p.is_update_way_R2L);
if(isreal(M))
[lambda, core_obj] = mexPalm4MSAReal(mex_params);
else
[lambda, core_obj] = mexPalm4MSACplx(mex_params);
end
F = Faust(core_obj, isreal(M));
end
%==========================================================================================
%> @brief Approximates M by A S_1 … S_n B using fact_palm4msa.
%>
%>
%> @Example
%> @code
%> import matfaust.fact.fact_palm4msa
%> import matfaust.factparams.*
%>
%> p = ParamsPalm4MSA(…
%> ConstraintList('spcol', 2, 10, 20, 'sp', 30, 20, 20),…
%> StoppingCriterion(50), 'is_verbose', false);
%> M = rand(10,10);
%> A = matfaust.rand(10,10);
%> B = matfaust.rand(20, 10);
%> [F, lamdba] = fact_palm4msa_constends(M, p, A, B)
%>
%> assert(norm(A - factors(F,1))/norm(A) <= eps(double(1)))
%> assert(norm(B - factors(F,4))/norm(B) <= eps(double(1)))
%>
%> @endcode
%==========================================================================================
function [F, lambda] = fact_palm4msa_constends(M, p, A, varargin)
import matfaust.factparams.*
import matfaust.fact.fact_palm4msa
if(~ ismatrix(A))
error('A must be a matrix.')
end
consA = ConstraintList('const', A, size(A,1), size(A,2));
new_consts = {};
new_consts = [ {consA.clist{:}}, {p.constraints{:}} ];
if(length(varargin) > 0)
B = varargin{1};
if(~ ismatrix(B))
error('B must be a matrix.')
end
consB = ConstraintList('const', B, size(B,1), size(B,2));
new_consts = [ new_consts, {consB.clist{:}} ];
end
new_consts = ConstraintList(new_consts{:});
p = ParamsPalm4MSA(new_consts, p.stop_crit, 'is_update_way_R2L', p.is_update_way_R2L, ...
'init_lambda', p.init_lambda, 'step_size', p.step_size, 'constant_step_size', ...
p.constant_step_size, 'is_verbose', p.is_verbose);
[F, lambda ] = fact_palm4msa(M, p);
f1 = factors(F, 1);
f1 = f1 / lambda;
nF = cell(1, numfactors(F));
nF{1} = f1;
for i=2:numfactors(F)
nF{i} = factors(F, i);
end
nF{2} = nF{2}*lambda;
F = matfaust.Faust(nF);
end
%==========================================================================================
%> @brief Computes the FGFT for the Laplacian matrix Lap.
%>
%>
%> @b Usage
%>
%> &nbsp;&nbsp;&nbsp; @b fgft_givens(Lap, J) calls the non-parallel Givens algorithm.<br/>
%> &nbsp;&nbsp;&nbsp; @b fgft_givens(Lap, J, 0) or fgft_givens(Lap, J, 1) do the same as in previous line.<br/>
%> &nbsp;&nbsp;&nbsp; @b fgft_givens(Lap, J, t) calls the parallel Givens algorithm (if t > 1, otherwise it calls basic Givens algorithm), see eigtj. <br/>
%> &nbsp;&nbsp;&nbsp; @b fgft_givens(Lap, J, t, 'verbosity', 2) same as above with a level of verbosity of 2 in output. <br/>
%>
%> @param Lap the Laplacian matrix as a numpy array. Must be real and symmetric.
%> @param J see eigtj
%> @param t see eigtj
%> @param verbosity see eigtj
%>
%> @retval [FGFT,D]:
%> - with FGFT being the Faust object representing the Fourier transform and,
%> - D as a sparse diagonal matrix of the eigenvalues in ascendant order along the rows/columns.
%>
%>
%> @b References:
%> - [1] Le Magoarou L., Gribonval R. and Tremblay N., "Approximate fast
%> graph Fourier transforms via multi-layer sparse approximations",
%> IEEE Transactions on Signal and Information Processing
%> over Networks 2018, 4(2), pp 407-420
%>
%>
%> <p> @b See @b also fact.eigtj, fact.fgft_palm
%>
%==========================================================================================
function [FGFT,D] = fgft_givens(Lap, J, varargin)
import matfaust.Faust
t = 1; % default value
verbosity = 0; % default value
if(~ ismatrix(Lap) || ~ isreal(Lap))
error('Lap must be a real matrix.')
end
if(size(Lap,1) ~= size(Lap,2))
error('Lap must be square')
end
if(~ isnumeric(J) || J-floor(J) > 0 || J <= 0)
error('J must be a positive integer.')
end
bad_arg_err = 'bad number of arguments.';
if(length(varargin) >= 1)
t = varargin{1};
if(~ isnumeric(t))
error('t must be a positive or nul integer.')
end
t = floor(abs(t));
t = min(t, J);
if(length(varargin) >= 2)
if(~ strcmp(varargin{2}, 'verbosity'))
error('arg. 4, if used, must be the str `verbosity''.')
end
if(length(varargin) == 3)
if(isnumeric(varargin{3}))
verbosity = floor(real(varargin{3}));
else
error('verbosity must be numeric')
end
else
error(bad_arg_err)
end
end
end
[core_obj, D] = mexfgftgivensReal(Lap, J, t, verbosity);
D = sparse(diag(D));
FGFT = Faust(core_obj, true);
end
%===================================================================================
%> @brief Computes the FGFT for the Fourier matrix U which should be the eigenvectors of the Laplacian Lap.
%>
%> @note this algorithm is a variant of fact_hierarchical.
%>
%> @param Lap The laplacian matrix.
%> @param U The Fourier matrix.
%> @param p The PALM hierarchical algorithm parameters.
%> @param init_D The initial diagonal vector. If none it will be the ones() vector by default.
%>
%> @retval [Uhat, Dhat, lambda, p]
%> - Uhat: the Faust factorization of U.
%> - Dhat: the diagonal matrix approximation of eigenvaules.
%> - lambda: see fact_hierarchical
%> - p: see fact_hierarchical
%>
%>
%> @b Example
%> @code
%> import matfaust.*
%> import matfaust.factparams.*
%> import matfaust.fact.fgft_palm
%>
%> % get the Laplacian
%> load('Laplacian_128_ring.mat');
%>
%> [U, D] = eig(Lap);
%> [D, I] = sort(diag(D));
%> D = diag(D);
%> U = U(:,I);
%>
%> dim = size(Lap, 1);
%>
%> nfacts = round(log2(dim)-3);
%> over_sp = 1.5; % sparsity overhead
%> dec_fact = .5; % decrease of the residuum sparsity
%>
%> % define the sparsity constraints for the factors
%> fact_cons = {};
%> res_cons = {};
%> for i=1:nfacts
%> fact_cons = [ fact_cons {ConstraintInt('sp', dim, dim, min(round(dec_fact^j*dim^2*over_sp), size(Lap,1)))} ];
%> res_cons = [ res_cons {ConstraintInt('sp', dim, dim, min(round(2*dim*over_sp), size(Lap, 1)))} ];
%> end
%>
%> % set the parameters for the PALM hierarchical algo.
%> params = ParamsHierarchicalFact(fact_cons, res_cons, StoppingCriterion(50), StoppingCriterion(100), 'step_size', 1e-6, 'constant_step_size', true, 'init_lambda', 1.0, 'is_fact_side_left', false);
%> %% compute FGFT for Lap, U, D
%> init_D_diag = diag(D);
%> [Uhat, Dhat, lambda, ~ ] = fgft_palm(U, Lap, params, init_D_diag);
%>
%> %% errors on FGFT and Laplacian reconstruction
%> err_U = norm(Uhat-U, 'fro')/norm(U, 'fro')
%> err_Lap = norm(Uhat*full(Dhat)*Uhat'-Lap, 'fro') / norm(Lap, 'fro')
%> % Output:
%> % Faust::HierarchicalFact<FPP,DEVICE,FPP2>::compute_facts : factorization 1/4
%> % Faust::HierarchicalFact<FPP,DEVICE,FPP2>::compute_facts : factorization 2/4
%> % Faust::HierarchicalFact<FPP,DEVICE,FPP2>::compute_facts : factorization 3/4
%> % Faust::HierarchicalFact<FPP,DEVICE,FPP2>::compute_facts : factorization 4/4
%> % err_U =
%> % 1.0013
%> % err_Lap =
%> % 0.9707
%>
%> @endcode
%>
%> <p> @b See @b also fact.fact_hierarchical, fact.eigtj, fact.fgft_givens
%>
%> @b References:
%> - [1] Le Magoarou L., Gribonval R. and Tremblay N., "Approximate fast
%> graph Fourier transforms via multi-layer sparse approximations",
%> IEEE Transactions on Signal and Information Processing
%> over Networks 2018, 4(2), pp 407-420
%> <https://hal.inria.fr/hal-01416110>
%> - [2] Le Magoarou L. and Gribonval R., "Are there approximate Fast
%> Fourier Transforms on graphs ?", ICASSP, 2016. <https://hal.inria.fr/hal-01254108>
%>
%===================================================================================
function varargout = fgft_palm(U, Lap, p, varargin)
import matfaust.Faust
import matfaust.factparams.*
import matfaust.fact.check_fact_mat
if(~ ismatrix(U) || ~ isnumeric(U) || ~ ismatrix(Lap) || ~ isnumeric(Lap))
error('U and Lap must be real or complex matrices.')
elseif(any(size(U) ~= size(Lap)) || any(size(Lap,1) ~= size(Lap,2)))
error('U and Lap must be square matrices of same size.')
end
% TODO: refactor with fact_hierarchical
if(length(varargin) == 1)
init_D = varargin{1};
if(~ ismatrix(init_D) || ~ isnumeric(init_D))
error('fgft_palm arg. 4 (init_D) must be a matrix')
end
if(size(init_D,1) ~= size(U,1))
error('fgft_palm arg. 4 (init_D) must be a diagonal vector of size == size(U,1).')
end
elseif(length(varargin) > 1)
error('fgft_palm, too many arguments.')
else % nargin == 0
init_D = ones(size(U,1),1);
if(~ isreal(U))
init_D = complex(init_D);
end
end
check_fact_mat('fgft_palm', U)
if(~ isa(p, 'ParamsHierarchicalFact') && ParamsFactFactory.is_a_valid_simplification(p))
p = ParamsFactFactory.createParams(U, p);
end
mex_constraints = cell(2, p.num_facts-1);
if(~ isa(p ,'ParamsHierarchicalFact'))
error('p must be a ParamsHierarchicalFact object.')
end
%mex_fact_constraints = cell(1, p.num_facts-1)
for i=1:p.num_facts-1
cur_cell = cell(1, 4);
cur_cell{1} = p.constraints{i}.name.conv2str();
cur_cell{2} = p.constraints{i}.param;
cur_cell{3} = p.constraints{i}.num_rows;
cur_cell{4} = p.constraints{i}.num_cols;
%mex_fact_constraints{i} = cur_cell;
mex_constraints{1,i} = cur_cell;
end
%mex_residuum_constraints = cell(1, p.num_facts-1)
for i=1:p.num_facts-1
cur_cell = cell(1, 4);
cur_cell{1} = p.constraints{i+p.num_facts-1}.name.conv2str();
cur_cell{2} = p.constraints{i+p.num_facts-1}.param;
cur_cell{3} = p.constraints{i+p.num_facts-1}.num_rows;
cur_cell{4} = p.constraints{i+p.num_facts-1}.num_cols;
%mex_residuum_constraints{i} = cur_cell;
mex_constraints{2,i} = cur_cell;
end
if(~ p.is_mat_consistent(U))
error('U''s number of columns must be consistent with the last residuum constraint defined in p. Likewise its number of rows must be consistent with the first factor constraint defined in p.')
end
% the setters for num_rows/cols verifies consistency with constraints
mex_params = struct('nfacts', p.num_facts, 'cons', {mex_constraints}, 'niter1', p.stop_crits{1}.num_its,'niter2', p.stop_crits{2}.num_its, 'sc_is_criterion_error', p.stop_crits{1}.is_criterion_error, 'sc_error_treshold', p.stop_crits{1}.error_treshold, 'sc_max_num_its', p.stop_crits{1}.max_num_its, 'sc_is_criterion_error2', p.stop_crits{2}.is_criterion_error, 'sc_error_treshold2', p.stop_crits{2}.error_treshold, 'sc_max_num_its2', p.stop_crits{2}.max_num_its, 'nrow', p.data_num_rows, 'ncol', p.data_num_cols, 'fact_side', p.is_fact_side_left, 'update_way', p.is_update_way_R2L, 'init_D', init_D, 'verbose', p.is_verbose, 'init_lambda', p.init_lambda);
if(isreal(U))
[lambda, core_obj, Ddiag] = mexHierarchical_factReal(U, mex_params, Lap);
else
[lambda, core_obj, Ddiag] = mexHierarchical_factCplx(U, mex_params, Lap);
end
D = sparse(diag(Ddiag));
F = Faust(core_obj, isreal(U));
varargout = {F, D, lambda, p};
end
%====================================================================
%> @brief Performs a singular value decomposition and returns the left and
%> right singular vectors as Faust transforms.
%>
%> @note this function is based on pyfaust.fact.eigtj.
%>
%> @param M: a real matrix.
%> @param J: see eigtj
%> @param t: see eigtj
%>
%> @retval [U,S,V]: U*full(S)*V' being the approximation of M.
%> - U: (sparse diagonal matrix) S the singular values in
%> descendant order.
%> - S: (Faust object) U the left-singular transform.
%> - V: (Faust object) V the right-singular transform.
%>
%> @Example
%> @code
%> % in a matlab terminal
%> >> import matfaust.fact.svdtj
%> >> M = rand(128,128)
%> >> [U,S,V] = svdtj(M,1024,64)
%> @endcode
%>
%====================================================================
function [U,S,V] = svdtj(M, J, varargin)
[W1,D1] = matfaust.FaustFactory.eigtj(M*M', J, varargin{:});
[W2,D2] = matfaust.FaustFactory.eigtj(M'*M, J, varargin{:});
S = diag(W1'*M*W2);
[~,I] = sort(abs(S), 'descend');
S = sparse(diag(S(I)));
sign_S = sign(S);
S = S*sign_S;
Id = eye(size(S));
U = W1(:,1:size(Id,1))*matfaust.Faust({Id(:,I),sign_S});
V = W2(:,1:size(Id,1))*matfaust.Faust(Id(:,I));
end
......@@ -32,63 +32,6 @@ classdef FaustFactory
end
methods(Static)
%==========================================================================================
%> @brief Factorizes the matrix M with Palm4MSA algorithm using the parameters set in p.
%>
%>
%> @param M the dense matrix to factorize.
%> @param p the ParamsPalm4MSA instance to define the algorithm parameters.
%>
%> @retval F the Faust object result of the factorization.
%> @retval [F, lambda] = fact_palm4msa(M, p) to optionally get lambda (scale).
%>
%> @b Example
%>
%> @code
%> import matfaust.*
%> import matfaust.factparams.*
%> M = rand(500, 32);
%> cons = ConstraintList('splin', 5, 500, 32, 'normcol', 1, 32, 32);
%> stop_crit = StoppingCriterion(200);
%> params = ParamsPalm4MSA(cons, stop_crit, 'is_update_way_R2L', false, 'init_lambda', 1.0);
%> F = FaustFactory.fact_palm4msa(M, params)
%> @endcode
%>
%> F =
%>
%> Faust size 500x32, density 0.22025, nnz_sum 3524, 2 factor(s):
%> - FACTOR 0 (real) SPARSE, size 500x32, density 0.15625, nnz 2500
%> - FACTOR 1 (real) SPARSE, size 32x32, density 1, nnz 1024
%>
%>
%==========================================================================================
function [F,lambda] = fact_palm4msa(M, p)
import matfaust.Faust
mex_constraints = cell(1, length(p.constraints));
matfaust.FaustFactory.check_fact_mat('FaustFactory.fact_palm4msa', M)
if(~ isa(p ,'matfaust.factparams.ParamsPalm4MSA'))
error('p must be a ParamsPalm4MSA object.')
end
for i=1:length(p.constraints)
cur_cell = cell(1, 4);
cur_cell{1} = p.constraints{i}.name.conv2str();
cur_cell{2} = p.constraints{i}.param;
cur_cell{3} = p.constraints{i}.num_rows;
cur_cell{4} = p.constraints{i}.num_cols;
mex_constraints{i} = cur_cell;
end
if(~ p.is_mat_consistent(M))
error('M''s number of columns must be consistent with the last residuum constraint defined in p. Likewise its number of rows must be consistent with the first factor constraint defined in p.')
end
% put mex_constraints in a cell array again because mex eats one level of array
mex_params = struct('data', M, 'nfacts', p.num_facts, 'cons', {mex_constraints}, 'init_facts', {p.init_facts}, 'niter', p.stop_crit.num_its, 'sc_is_criterion_error', p.stop_crit.is_criterion_error, 'sc_error_treshold', p.stop_crit.error_treshold, 'sc_max_num_its', p.stop_crit.max_num_its, 'update_way', p.is_update_way_R2L);
if(isreal(M))
[lambda, core_obj] = mexPalm4MSAReal(mex_params);
else
[lambda, core_obj] = mexPalm4MSACplx(mex_params);
end
F = Faust(core_obj, isreal(M));
end
%==========================================================================================
......@@ -379,312 +322,7 @@ classdef FaustFactory
varargout = {F, lambda, p};
end
%===================================================================================
%> @brief Computes the FGFT for the Fourier matrix U which should be the eigenvectors of the Laplacian Lap.
%>
%> @note this algorithm is a variant of FaustFactory.fact_hierarchical.
%>
%> @param Lap The laplacian matrix.
%> @param U The Fourier matrix.
%> @param p The PALM hierarchical algorithm parameters.
%> @param init_D The initial diagonal vector. If none it will be the ones() vector by default.
%>
%> @retval [Uhat, Dhat, lambda, p]
%> - Uhat: the Faust factorization of U.
%> - Dhat: the diagonal matrix approximation of eigenvaules.
%> - lambda: see FaustFactory.fact_hierarchical
%> - p: see FaustFactory.fact_hierarchical
%>
%>
%> @b Example
%> @code
%> import matfaust.*
%> import matfaust.factparams.*
%>
%> % get the Laplacian
%> load('Laplacian_128_ring.mat');
%>
%> [U, D] = eig(Lap);
%> [D, I] = sort(diag(D));
%> D = diag(D);
%> U = U(:,I);
%>
%> dim = size(Lap, 1);
%>
%> nfacts = round(log2(dim)-3);
%> over_sp = 1.5; % sparsity overhead
%> dec_fact = .5; % decrease of the residuum sparsity
%>
%> % define the sparsity constraints for the factors
%> fact_cons = {};
%> res_cons = {};
%> for i=1:nfacts
%> fact_cons = [ fact_cons {ConstraintInt('sp', dim, dim, min(round(dec_fact^j*dim^2*over_sp), size(Lap,1)))} ];
%> res_cons = [ res_cons {ConstraintInt('sp', dim, dim, min(round(2*dim*over_sp), size(Lap, 1)))} ];
%> end
%>
%> % set the parameters for the PALM hierarchical algo.
%> params = ParamsHierarchicalFact(fact_cons, res_cons, StoppingCriterion(50), StoppingCriterion(100), 'step_size', 1e-6, 'constant_step_size', true, 'init_lambda', 1.0, 'is_fact_side_left', false);
%> %% compute FGFT for Lap, U, D
%> init_D_diag = diag(D);
%> [Uhat, Dhat, lambda, ~ ] = FaustFactory.fgft_palm(U, Lap, params, init_D_diag);
%>
%> %% errors on FGFT and Laplacian reconstruction
%> err_U = norm(Uhat-U, 'fro')/norm(U, 'fro')
%> err_Lap = norm(Uhat*full(Dhat)*Uhat'-Lap, 'fro') / norm(Lap, 'fro')
%> % Output:
%> % Faust::HierarchicalFact<FPP,DEVICE,FPP2>::compute_facts : factorization 1/4
%> % Faust::HierarchicalFact<FPP,DEVICE,FPP2>::compute_facts : factorization 2/4
%> % Faust::HierarchicalFact<FPP,DEVICE,FPP2>::compute_facts : factorization 3/4
%> % Faust::HierarchicalFact<FPP,DEVICE,FPP2>::compute_facts : factorization 4/4
%> % err_U =
%> % 1.0013
%> % err_Lap =
%> % 0.9707
%>
%> @endcode
%>
%> <p> @b See @b also FaustFactory.fact_hierarchical, FaustFactory.eigtj, FaustFactory.fgft_givens
%>
%> @b References:
%> - [1] Le Magoarou L., Gribonval R. and Tremblay N., "Approximate fast
%> graph Fourier transforms via multi-layer sparse approximations",
%> IEEE Transactions on Signal and Information Processing
%> over Networks 2018, 4(2), pp 407-420
%> <https://hal.inria.fr/hal-01416110>
%> - [2] Le Magoarou L. and Gribonval R., "Are there approximate Fast
%> Fourier Transforms on graphs ?", ICASSP, 2016. <https://hal.inria.fr/hal-01254108>
%>
%===================================================================================
function varargout = fgft_palm(U, Lap, p, varargin)
import matfaust.Faust
import matfaust.factparams.*
if(~ ismatrix(U) || ~ isnumeric(U) || ~ ismatrix(Lap) || ~ isnumeric(Lap))
error('U and Lap must be real or complex matrices.')
elseif(any(size(U) ~= size(Lap)) || any(size(Lap,1) ~= size(Lap,2)))
error('U and Lap must be square matrices of same size.')
end
% TODO: refactor with fact_hierarchical
if(length(varargin) == 1)
init_D = varargin{1};
if(~ ismatrix(init_D) || ~ isnumeric(init_D))
error('fgft_palm arg. 4 (init_D) must be a matrix')
end
if(size(init_D,1) ~= size(U,1))
error('fgft_palm arg. 4 (init_D) must be a diagonal vector of size == size(U,1).')
end
elseif(length(varargin) > 1)
error('fgft_palm, too many arguments.')
else % nargin == 0
init_D = ones(size(U,1),1);
if(~ isreal(U))
init_D = complex(init_D);
end
end
matfaust.FaustFactory.check_fact_mat('FaustFactory.fgft_palm', U)
if(~ isa(p, 'ParamsHierarchicalFact') && ParamsFactFactory.is_a_valid_simplification(p))
p = ParamsFactFactory.createParams(U, p);
end
mex_constraints = cell(2, p.num_facts-1);
if(~ isa(p ,'ParamsHierarchicalFact'))
error('p must be a ParamsHierarchicalFact object.')
end
%mex_fact_constraints = cell(1, p.num_facts-1)
for i=1:p.num_facts-1
cur_cell = cell(1, 4);
cur_cell{1} = p.constraints{i}.name.conv2str();
cur_cell{2} = p.constraints{i}.param;
cur_cell{3} = p.constraints{i}.num_rows;
cur_cell{4} = p.constraints{i}.num_cols;
%mex_fact_constraints{i} = cur_cell;
mex_constraints{1,i} = cur_cell;
end
%mex_residuum_constraints = cell(1, p.num_facts-1)
for i=1:p.num_facts-1
cur_cell = cell(1, 4);
cur_cell{1} = p.constraints{i+p.num_facts-1}.name.conv2str();
cur_cell{2} = p.constraints{i+p.num_facts-1}.param;
cur_cell{3} = p.constraints{i+p.num_facts-1}.num_rows;
cur_cell{4} = p.constraints{i+p.num_facts-1}.num_cols;
%mex_residuum_constraints{i} = cur_cell;
mex_constraints{2,i} = cur_cell;
end
if(~ p.is_mat_consistent(U))
error('U''s number of columns must be consistent with the last residuum constraint defined in p. Likewise its number of rows must be consistent with the first factor constraint defined in p.')
end
% the setters for num_rows/cols verifies consistency with constraints
mex_params = struct('nfacts', p.num_facts, 'cons', {mex_constraints}, 'niter1', p.stop_crits{1}.num_its,'niter2', p.stop_crits{2}.num_its, 'sc_is_criterion_error', p.stop_crits{1}.is_criterion_error, 'sc_error_treshold', p.stop_crits{1}.error_treshold, 'sc_max_num_its', p.stop_crits{1}.max_num_its, 'sc_is_criterion_error2', p.stop_crits{2}.is_criterion_error, 'sc_error_treshold2', p.stop_crits{2}.error_treshold, 'sc_max_num_its2', p.stop_crits{2}.max_num_its, 'nrow', p.data_num_rows, 'ncol', p.data_num_cols, 'fact_side', p.is_fact_side_left, 'update_way', p.is_update_way_R2L, 'init_D', init_D, 'verbose', p.is_verbose, 'init_lambda', p.init_lambda);
if(isreal(U))
[lambda, core_obj, Ddiag] = mexHierarchical_factReal(U, mex_params, Lap);
else
[lambda, core_obj, Ddiag] = mexHierarchical_factCplx(U, mex_params, Lap);
end
D = sparse(diag(Ddiag));
F = Faust(core_obj, isreal(U));
varargout = {F, D, lambda, p};
end
%==========================================================================================
%> @brief Computes the FGFT for the Laplacian matrix Lap.
%>
%>
%> @b Usage
%>
%> &nbsp;&nbsp;&nbsp; @b fgft_givens(Lap, J) calls the non-parallel Givens algorithm.<br/>
%> &nbsp;&nbsp;&nbsp; @b fgft_givens(Lap, J, 0) or fgft_givens(Lap, J, 1) do the same as in previous line.<br/>
%> &nbsp;&nbsp;&nbsp; @b fgft_givens(Lap, J, t) calls the parallel Givens algorithm (if t > 1, otherwise it calls basic Givens algorithm), see FaustFactory.eigtj. <br/>
%> &nbsp;&nbsp;&nbsp; @b fgft_givens(Lap, J, t, 'verbosity', 2) same as above with a level of verbosity of 2 in output. <br/>
%>
%> @param Lap the Laplacian matrix as a numpy array. Must be real and symmetric.
%> @param J see FaustFactory.eigtj
%> @param t see FaustFactory.eigtj
%> @param verbosity see FaustFactory.eigtj
%>
%> @retval [FGFT,D]:
%> - with FGFT being the Faust object representing the Fourier transform and,
%> - D as a sparse diagonal matrix of the eigenvalues in ascendant order along the rows/columns.
%>
%>
%> @b References:
%> - [1] Le Magoarou L., Gribonval R. and Tremblay N., "Approximate fast
%> graph Fourier transforms via multi-layer sparse approximations",
%> IEEE Transactions on Signal and Information Processing
%> over Networks 2018, 4(2), pp 407-420
%>
%>
%> <p> @b See @b also FaustFactory.eigtj, FaustFactory.fgft_palm
%>
%==========================================================================================
function [FGFT,D] = fgft_givens(Lap, J, varargin)
import matfaust.Faust
t = 1; % default value
verbosity = 0; % default value
if(~ ismatrix(Lap) || ~ isreal(Lap))
error('Lap must be a real matrix.')
end
if(size(Lap,1) ~= size(Lap,2))
error('Lap must be square')
end
if(~ isnumeric(J) || J-floor(J) > 0 || J <= 0)
error('J must be a positive integer.')
end
bad_arg_err = 'bad number of arguments.';
if(length(varargin) >= 1)
t = varargin{1};
if(~ isnumeric(t))
error('t must be a positive or nul integer.')
end
t = floor(abs(t));
t = min(t, J);
if(length(varargin) >= 2)
if(~ strcmp(varargin{2}, 'verbosity'))
error('arg. 4, if used, must be the str `verbosity''.')
end
if(length(varargin) == 3)
if(isnumeric(varargin{3}))
verbosity = floor(real(varargin{3}));
else
error('verbosity must be numeric')
end
else
error(bad_arg_err)
end
end
end
[core_obj, D] = mexfgftgivensReal(Lap, J, t, verbosity);
D = sparse(diag(D));
FGFT = Faust(core_obj, true);
end
%==========================================================================================
%> @brief Computes the eigenvalues and the eigenvectors transform (as a Faust object) using the truncated Jacobi algorithm.
%>
%> The eigenvalues and the eigenvectors are approximate. The trade-off between accuracy and sparsity can be set through the parameters J and t.
%>
%>
%> @b Usage
%>
%> &nbsp;&nbsp;&nbsp; @b eigtj(M, J) calls the non-parallel Givens algorithm.<br/>
%> &nbsp;&nbsp;&nbsp; @b eigtj(M, J, 0) or eigtj(M, J, 1) do the same as in previous line.<br/>
%> &nbsp;&nbsp;&nbsp; @b eigtj(M, J, t) calls the parallel Givens algorithm (if t > 1, otherwise it calls basic Givens algorithm)<br/>
%> &nbsp;&nbsp;&nbsp; @b eigtj(M, J, t, 'verbosity', 2) same as above with a level of verbosity of 2 in output. <br/>
%>
%> @param M the matrix to diagonalize. Must be real and symmetric.
%> @param J 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).
%> @param t 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 > size(M,1)/2 won't lead to the desired effect because the maximum number of rotation matrices per factor is anyway size(M,1)/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.
%> @param verbosity the level of verbosity, the greater the value the more info. is displayed.
%>
%> @retval [V,D]
%> - 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 eigenvalues set in D.
%> - D the approximate sparse diagonal matrix of the eigenvalues (in ascendant order along the rows/columns).
%>
%> @b Example
%> @code
%> import matfaust.*
%>
%> % get a Laplacian to diagonalize
%> load('Laplacian_256_community.mat')
%> % do it
%> [Uhat, Dhat] = FaustFactory.eigtj(Lap, size(Lap,1)*100, size(Lap, 1)/2, 'verbosity', 2)
%> % Uhat is the Fourier matrix/eigenvectors approximattion as a Faust (200 factors + permutation mat.)
%> % Dhat the eigenvalues diagonal matrix approx.
%> @endcode
%>
%>
%>
%> @b References:
%> - [1] Le Magoarou L., Gribonval R. and Tremblay N., "Approximate fast
%> graph Fourier transforms via multi-layer sparse approximations",
%> IEEE Transactions on Signal and Information Processing
%> over Networks 2018, 4(2), pp 407-420
%>
%> <p> @b See @b also FaustFactory.fgft_givens, FaustFactory.fgft_palm
%>
%==========================================================================================
function [V,D] = eigtj(M, J, varargin)
[V, D] = matfaust.FaustFactory.fgft_givens(M, J, varargin{:});
V = factors(V,1:numfactors(V))
% copy seems unecessary but it's to workaround a bug (temporarily)
end
%====================================================================
%> @brief Performs a singular value decomposition and returns the left and
%> right singular vectors as Faust transforms.
%>
%> @note this function is based on FaustFactory.eigtj.
%>
%> @param M: a real matrix.
%> @param J: see FaustFactory.eigtj
%> @param t: see FaustFactory.eigtj
%>
%> @retval [U,S,V]: U*full(S)*V' being the approximation of M.
%> - U: (sparse diagonal matrix) S the singular values in
%> descendant order.
%> - S: (Faust object) U the left-singular transform.
%> - V: (Faust object) V the right-singular transform.
%>
%> @Example
%> @code
%> % in a matlab terminal
%> >> import matfaust.*
%> >> M = rand(128,128)
%> >> [U,S,V] = FaustFactory.svdtj(M,1024,64)
%> @endcode
%>
%====================================================================
function [U,S,V] = svdtj(M, J, varargin)
[W1,D1] = matfaust.FaustFactory.eigtj(M*M', J, varargin{:});
[W2,D2] = matfaust.FaustFactory.eigtj(M'*M, J, varargin{:});
S = diag(W1'*M*W2);
[~,I] = sort(abs(S), 'descend');
S = sparse(diag(S(I)));
sign_S = sign(S);
S = S*sign_S;
Id = eye(size(S));
U = W1(:,1:size(Id,1))*matfaust.Faust({Id(:,I),sign_S});
V = W2(:,1:size(Id,1))*matfaust.Faust(Id(:,I));
end
end
methods(Access = private, Static)
......
This diff is collapsed.
This diff is collapsed.
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment