-
hhakim authored
Add FGFT error benchmark scripts (python and matlab) and a graph generation script using GSP matlab toolbox. The goal of these scripts is to validate the C++ port of the FGFT algorithms and their python wrapper against the reference Matlab prototypes from FAuST 1.03. The two FGFT algorithms are tested : PALM FGFT and Givens FGFT.
hhakim authoredAdd FGFT error benchmark scripts (python and matlab) and a graph generation script using GSP matlab toolbox. The goal of these scripts is to validate the C++ port of the FGFT algorithms and their python wrapper against the reference Matlab prototypes from FAuST 1.03. The two FGFT algorithms are tested : PALM FGFT and Givens FGFT.
benchmark_Lap_diag.py 8.57 KiB
from scipy.io import loadmat
from sys import argv
from os.path import dirname, sep, join, exists
from numpy import diag, copy, log2, count_nonzero, size, loadtxt, savetxt
import numpy, re
from numpy.linalg import norm
from pyfaust.factparams import *
from pyfaust import *
from pylab import *
if __name__ == '__main__':
data_file = join(dirname(argv[0]),
'../../../data/mat/Laps_U_Ds-6_graph_types-dims_128_to_1024.mat')
if(not exists(data_file)):
raise Exception(data_file+" not found. You need to run"
" misc/test/src/Matlab/get_graph_test_UDLap_matrices.m"
" for generation.")
matfile = loadmat(data_file)
Laps = matfile['Laps'][0]
Ds = matfile['Ds'][0]
Us = matfile['Us'][0]
start_i=0
num_laps = Laps.shape[0]
num_laps = 45
hier_fgft_errs = empty(num_laps)
givens_errs = empty(num_laps)
par_givens_errs = empty(num_laps)
U_hier_fgft_errs = empty(num_laps)
U_givens_errs = empty(num_laps)
U_par_givens_errs = empty(num_laps)
if exists('benchmark_pyfaust_fgft.txt'):
saved_lap_errs = loadtxt('benchmark_pyfaust_fgft.txt')
len_saved = saved_lap_errs[0].shape[0]
givens_errs[:len_saved], par_givens_errs[:len_saved],\
hier_fgft_errs[:len_saved] = saved_lap_errs
saved_U_errs = loadtxt('benchmark_pyfaust_fgft_U.txt')
U_givens_errs[:len_saved], U_par_givens_errs[:len_saved], \
U_hier_fgft_errs[:len_saved] = \
saved_U_errs
start_i = len_saved
if(saved_U_errs[0].shape[0] != saved_lap_errs[0].shape[0]):
raise Exception("Error files not consistent, delete and recompute")
for i in range(start_i, num_laps):
print('Lap ',i,'/',num_laps)
U,D,Lap = Us[i],Ds[i],Laps[i].astype(numpy.float)
dim = Lap.shape[0]
nfacts = int(round(log2(dim))-3)
over_sp = 1.5 # sparsity overhead
dec_fact = .5 # decrease of the residum sparsity
fact_cons, res_cons = [], []
for j in range(1, nfacts):
fact_cons += [ ConstraintInt('sp', dim, dim,
min(int(round(dec_fact**j*dim**2*over_sp)),size(Lap)))
]
res_cons += [ ConstraintInt('sp', dim, dim,
min(int(round(2*dim*over_sp)),size(Lap))) ]
params = ParamsHierarchicalFact(fact_cons, res_cons,
StoppingCriterion(num_its=50),
StoppingCriterion(num_its=100),
step_size=1.0000e-06,
constant_step_size=False,
init_lambda=1.0,
is_fact_side_left=False)
F = FaustFactory.fact_hierarchical(U,params)
print(F)
complexity_global = F.nnz_sum()
rc_palm = complexity_global / count_nonzero(U)
diag_init_D = diag(D)
diag_init_D = numpy.copy(diag_init_D)
F, Dhat = FaustFactory.fgft_palm(U, Lap,
params,
diag_init_D)
hier_fgft_err = norm((F.todense()*diag(Dhat))*F.T.todense()-Lap,"fro")/norm(Lap,"fro")
hier_fgft_errs[i] = hier_fgft_err
U_hier_fgft_errs[i] = (F-U).norm("fro")/norm(U,"fro")
#nfacts = complexity_global/(2*dim)
#J = round(nfacts*(dim/2))
J = round(complexity_global/4)
F, Dhat = FaustFactory.fgft_givens(Lap, J, 0)
givens_err = norm((F*Dhat.todense())*F.T.todense()-Lap,'fro')/norm(Lap,'fro')
givens_errs[i] = givens_err
U_givens_errs[i] = (F-U).norm("fro")/norm(U,"fro")
F, Dhat = FaustFactory.fgft_givens(Lap, J, int(dim/2))
print("J=", J)
print("nnz_sum FGFT givens parallel=", F.nnz_sum(), "num of facts:",
F.get_num_factors())
par_givens_err = norm((F*Dhat.todense())*F.T.todense()-Lap,'fro')/norm(Lap,'fro')
par_givens_errs[i] = par_givens_err
U_par_givens_errs[i] = (F-U).norm("fro")/norm(U,"fro")
print("err:", norm(Lap * F.toarray() - F.toarray()*Dhat)/norm(Lap))
savetxt('benchmark_pyfaust_fgft.txt', np.array([givens_errs[:i+1],
par_givens_errs[:i+1],
hier_fgft_errs[:i+1]]))
savetxt('benchmark_pyfaust_fgft_U.txt', np.array([U_givens_errs[:i+1],
U_par_givens_errs[:i+1],
U_hier_fgft_errs[:i+1]]))
matlab_matfile = None
if(len(argv) > 1 and exists(argv[1]) and re.match(".*.mat", argv[1])):
matlab_matfile = loadmat(argv[1])
matlab_hier_fgft_errs = [matlab_matfile['hier_errs'][0,i][0,0] for i in
range(matlab_matfile['hier_errs'].shape[1]) ]
matlab_givens_errs = [matlab_matfile['givens_errs'][0,i][0,0] for i in
range(matlab_matfile['givens_errs'].shape[1]) ]
matlab_par_givens_errs = [matlab_matfile['par_givens_errs'][0,i][0,0] for i in
range(matlab_matfile['par_givens_errs'].shape[1]) ]
print('hier_fgft_errs = ', hier_fgft_errs)
print('givens_errs = ', givens_errs)
print('par_givens_errs = ', par_givens_errs)
print('U_hier_fgft_errs = ', U_hier_fgft_errs)
print('U_givens_errs = ', U_givens_errs)
print('U_par_givens_errs = ', U_par_givens_errs)
plt.rcParams['figure.figsize'] = [18.0, 12]
lap_indices = arange(num_laps)
grid(True, which='both', axis='both')
xticks(lap_indices)
scatter(lap_indices, givens_errs, c='b', marker='o', label='pyfaust_givens_errs')
scatter(lap_indices, par_givens_errs, c='r', marker='o', label='pyfaust_par_givens_errs')
scatter(lap_indices, hier_fgft_errs, c='y', marker='o', label='pyfaust_hier_fgft_errs')
if matlab_matfile:
scatter(lap_indices, matlab_hier_fgft_errs[:len(lap_indices)],
label='matlab_hier_fgft_errs', c='y', marker='+', s=315)
scatter(lap_indices, matlab_givens_errs[:len(lap_indices)],
label='matlab_givens_errs', c='b', marker='+', s=315)
scatter(lap_indices, matlab_par_givens_errs[:len(lap_indices)],
label='matlab_par_givens_errs', c='r', marker='+', s=315)
ylabel("Relative error on Laplacian (reconstructed with FGFT)")
xlabel("Laplacians\n i-th Lap. is for"
" erdos-renyi if i==0 (mod 6)\ncommunity if i==1 (mod 6)\nsensor if i=="
"2 (mod 6)\npath if i==3 (mod 6)\nrandom ring if i==4 (mod6), ring if i=="
"5 (mod6)")
title("Error benchmark on "+str(num_laps)+" Laplacians (128x128)")
legend()
savefig('benchmark_Lap_diag_pyfaust_Laplacian_figure.png')
figure(2)
grid(True, which='both', axis='both')
if matlab_matfile:
matlab_U_hier_fgft_errs = [matlab_matfile['U_hier_errs'][0,i][0,0] for i in
range(matlab_matfile['U_hier_errs'].shape[1]) ]
matlab_U_givens_errs = [matlab_matfile['U_givens_errs'][0,i][0,0] for i in
range(matlab_matfile['U_givens_errs'].shape[1]) ]
matlab_U_par_givens_errs = [matlab_matfile['U_par_givens_errs'][0,i][0,0] for i in
range(matlab_matfile['U_par_givens_errs'].shape[1]) ]
scatter(lap_indices, matlab_U_hier_fgft_errs[:len(lap_indices)],
label='matlab_hier_fgft_errs', c='y', marker='+', s=315)
scatter(lap_indices, matlab_U_givens_errs[:len(lap_indices)],
label='matlab_givens_errs', c='b', marker='+', s=315)
scatter(lap_indices, matlab_U_par_givens_errs[:len(lap_indices)],
label='matlab_par_givens_errs', c='r', marker='+', s=315)
xticks(lap_indices)
scatter(lap_indices, U_givens_errs, c='b', marker='o', label='pyfaust_givens_errs')
scatter(lap_indices, U_par_givens_errs, c='r', marker='o', label='pyfaust_par_givens_errs')
scatter(lap_indices, U_hier_fgft_errs, c='y', marker='o', label='pyfaust_hier_fgft_errs')
ylabel("Relative error of FGFT versus U (Fourier matrix)")
xlabel("Fourier matrices")
title("Error benchmark on "+str(num_laps)+" Graph Fourier matrices (128x128)")
legend()
if matlab_matfile:
print('matlab_hier_fgft_errs = ', matlab_hier_fgft_errs)
print('matlab_givens_errs = ', matlab_givens_errs)
print('matlab_par_givens_errs = ', matlab_par_givens_errs)
print('matlab_U_hier_fgft_errs = ', matlab_U_hier_fgft_errs)
print('matlab_U_givens_errs = ', matlab_U_givens_errs)
print('matlab_U_par_givens_errs = ', matlab_U_par_givens_errs)
savefig('benchmark_Lap_diag_pyfaust_Fourier_figure.png')
show()