%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%% test/driver routine
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear all
load A.mat % 802 x 802 test matrix

n = size(A,1);
r = 16; % HSS block row size
tol = 1e-6;
[tr,m] = npart(n,r);  % binary tree; partition sequence of n

%%%%%%%%%% Generation of HSS approximation
[D,U,R,B,W,V,nflops1] = mat2hss(A,tr,m,'tol',tol);

%%%%%%%%%% HSS system solution test
x0 = randn(n,1); b = A*x0;
[x,nflops2] = hssulvsol(tr,D,U,R,B,W,V,length(tr),b);
fprintf('relative residual: %e\n',norm(A*x-b)/norm(b));
fprintf('HSS construction flops: %e\nHSS solution     flops: %e\n',nflops1,nflops2);
fprintf('Regular solution flops: %e\n',2/3*n^3);

%%%%%%%%%% Test the accuracy of the HSS matrix approximation
A1 = hss2mat(D,U,R,B,W,V,tr);
fprintf('||A-tilde A||_2/||A||_2: %e\n',norm(A-A1)/norm(A));