Commit 4692ae65 authored by Rafael Laboissiere's avatar Rafael Laboissiere

Imported Upstream version 1.0.1

parents
This diff is collapsed.
Name: sparsersb
Version: 1.0.1
Date: 2016-08-01
Author: Michele Martone <michelemartone@users.sourceforge.net>
Maintainer: Michele Martone <michelemartone@users.sourceforge.net>
Title: Interface to the librsb package implementing the RSB sparse matrix format.
Description: Interface to the librsb package implementing the RSB sparse matrix format for fast shared-memory sparse matrix computations.
Depends: octave (>= 4.0.0)
License: GPLv3+
Url: http://librsb.sourceforge.net/
Categories: Sparse Matrix Computations
BuildRequires: librsb
sparsersb >> Interface to the librsb package implementing the RSB sparse matrix format.
Support
sparsersb
Summary of important user-visible changes for releases of the sparsersb package
===============================================================================
sparsersb-1.0.1 Release Date: 2016-08-05
===============================================================================
** Thought to be used with the latest librsb-1.2.0.
** Changed sparsersb's `configure --help':
** - options to build librsb from a tarball (via configure or LIBRSB_TARBALL)
** - options to use librsb-config
** - options to override librsb-config
** - you can override the default C++11 flag
===============================================================================
sparsersb-1.0.0 Release Date: 2015-05-31
===============================================================================
** First Packaged Release. Intended to work with librsb-1.2.
#
# Copyright (C) 2011-2015 Michele Martone <michelemartone _AT_ users.sourceforge.net>
#
# This program is free software; you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation; either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program; if not, see <http://www.gnu.org/licenses/>.
#
# This documentation is yet to be written.
disp "This demo is to be written"
#
# Copyright (C) 2011-2015 Michele Martone <michelemartone _AT_ users.sourceforge.net>
#
# This program is free software; you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation; either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program; if not, see <http://www.gnu.org/licenses/>.
#
#
# Linear Solvers benchmarks using sparsersb.
#
# TODO: this file shall host some linear system solution benchmarks using sparsersb.
# It may serve as a reference point when profiling sparsersb/librsb.
# Please note that sparsersb is optimized for large matrices.
#
1; # This is a script
function lsb_compare(A)
n=rows(A);
maxit = n;
b = ones (n, 1);
P = diag (diag (A));
[i,j,v]=find(sparse(A));
minres=1e-7;
printf("Solving a system of %d equations, %d nonzeroes.\n",n,nnz(A));
tic; Ao = sparse (i,j,v,n,n);obt=toc;
onz=nnz(Ao);
tic; [X, FLAG, RELRES, ITER] = gmres (Ao, b, [], minres, maxit, P); odt=toc;
cs="Octave ";
onv=norm(Ao*X-b);
oRELRES=RELRES;
printf("%s took %.4f = %.4f + %.4f s and gave residual %g, flag %d, error norm %g.\n",cs,obt+odt,obt,odt,RELRES,FLAG,onv);
tic; Ar = sparsersb (i,j,v,n,n);rbt=toc;
#tic; Ar = sparsersb (Ao);rbt=toc;
rnz=nnz(Ar);
tic; [X, FLAG, RELRES, ITER] = gmres (Ar, b, [], minres, maxit, P); rdt=toc;
cs="sparsersb";
rnv=norm(Ar*X-b);
printf("%s took %.4f = %.4f + %.4f s and gave residual %g, flag %d, error norm %g.\n",cs,rbt+rdt,rbt,rdt,RELRES,FLAG,rnv);
if (onz != rnz)
printf("Error: seems like matrices don't match: %d vs %d nonzeroes!\n",onz,rnz);
quit(1);
else
end
if (RELRES>minres ) && (oRELRES<minres )
printf("Error: sparsersb was not able to solve a system octave did (residuals: %g vs %g)!",RELRES,oRELRES);
quit(1);
else
printf("Both systems were solved, speedups for overall: %g, constructor: %g, iterations: %g.\n",(obt+odt)/(rbt+rdt),(obt)/(rbt),(odt)/(rdt));
end
printf("\n");
end
# This one is based on what Carlo De Falco posted on the octave-dev mailing list:
# (he used n=1000, k=15)
n = 4;
k = 1;
A= k * eye (n) + sprandn (n, n, .8);
lsb_compare(A);
n = 100;
k = 5;
A= k * eye (n) + sprandn (n, n, .2);
lsb_compare(A);
n = 1000;
k = 15;
A= k * eye (n) + sprandn (n, n, .2);
lsb_compare(A);
n = 5000;
k = 1500;
A= k * eye (n) + sprandn (n, n, .2);
lsb_compare(A);
#nx=40,ny=20;
nx0=10;ny0=10;
for xm=1:2
#for ym=1:2
for ym=1:1
nx=nx0*(2^xm),ny=ny0*(2^ym);
hx=1/(nx+1);
hy=1/(ny+1);
# a symmetric example, from andreas stahel's notes on solving linear systems
Dxx=spdiags([ones(nx,1)-2*ones(nx,1) ones(nx,1)],[-1 0 1],nx,nx)/(hx^2);
Dyy=spdiags([ones(ny,1)-2*ones(ny,1) ones(ny,1)],[-1 0 1],ny,ny)/(hy^2);
A = - kron(Dxx, speye(ny))-kron(speye(nx),Dyy);
[xx,yy]=meshgrid(linspace(hx,1-hx,nx),linspace(ny,1-hy,ny));
b=xx(:);
lsb_compare(A);
# a symmetric example, from andreas stahel's notes on solving linear systems
Dy=spdiags([ones(ny,1) ones(ny,1)],[-1 1],ny,ny)/(2*hy);
A = - kron(Dxx, speye(ny))-kron(speye(nx),Dyy) + 5*kron(speye(nx),Dy);
[xx,yy]=meshgrid(linspace(hx,1-hx,nx),linspace(ny,1-hy,ny));
b=xx(:);
lsb_compare(A);
end
end
printf "All done."
quit(1);
#
# Copyright (C) 2011-2015 Michele Martone <michelemartone _AT_ users.sourceforge.net>
#
# This program is free software; you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation; either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program; if not, see <http://www.gnu.org/licenses/>.
#
1; # This is a script
# once complete, this program will benchmark our matrix implementation against octave's
cmt="#";
#for n_=1:6*0+1
for n_=1:6
for ro=0:1
n=n_*1000;
m=k=n;
# making vectors
b=linspace(1,1,n)';
ox=linspace(1,1,n)';
bx=linspace(1,1,n)';
# making matrices
r=(rand(n)>.6);
om=sparse(r);
nz=nnz(om);
M=10^6;
if ro==1
printf("%s%s\n",cmt,"reordering with colamd");
pct=-time;
p=colamd(om);
pct+=time;
pat=-time;
om=om(:,p);
pat+=time;
# TODO: use an array to select/specify the different reordering algorithms
printf("%g\t%g\t(%s)\n",(nz/M)/pct,(nz/M)/pat,"mflops for pct/pat");
else
printf("%s%s\n",cmt,"no reordering");
end
#bm=sparsevbr(om);
bm=sparsersb(sparse(om));
#bm=sparsersb3(sparse(om));
# stats
flops=2*nz;
## spmv
ot=-time; ox=om*b; ot+=time;
#
bt=-time; bx=bm*b; bt+=time;
t=ot; p=["octave-",version]; mflops=(flops/M)/t;
printf("%s\t%d\t%d\t%d\t%g\t%s\n","*",m,k,nz,mflops, p);
t=bt; p=["RSB"]; mflops=(flops/M)/t;
printf("%s\t%d\t%d\t%d\t%g\t%s\n","*",m,k,nz,mflops, p);
## spmvt
ot=-time; ox=om.'*b; ot+=time;
#
bt=-time; bx=bm.'*b; bt+=time;
t=ot; p=["octave-",version]; mflops=(flops/M)/t;
printf("%s\t%d\t%d\t%d\t%g\t%s\n","*",m,k,nz,mflops, p);
t=bt; p=["RSB"]; mflops=(flops/M)/t;
printf("%s\t%d\t%d\t%d\t%g\t%s\n","*",m,k,nz,mflops, p);
## spgemm
ot=-time; ox=om*om; ot+=time;
#
bt=-time; bx=bm*bm; bt+=time;
t=ot; p=["octave-",version]; mflops=n*(flops/M)/t;
printf("%s\t%d\t%d\t%d\t%g\t%s\n","OCT_SPGEMM",m,k,nz,mflops, p);
t=bt; p=["RSB"]; mflops=n*(flops/M)/t;
printf("%s\t%d\t%d\t%d\t%g\t%s\n","RSB_SPGEMM",m,k,nz,mflops, p);
endfor
endfor
#!/usr/bin/octave -q
#
# Copyright (C) 2011-2015 Michele Martone <michelemartone _AT_ users.sourceforge.net>
#
# This program is free software; you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation; either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program; if not, see <http://www.gnu.org/licenses/>.
#
1; # This is a script
# a benchmark program for octave/matlab
# TODO: fix output format
# TODO: correct symmetric / hermitian matrices handling
# TODO: sound, time-and-runs-based benchmarking criteria
n=10;
function printbenchline(matrixname,opname,sw,times,nnz,tottime,mxxops,bpnz,msstr)
printf("FIXME (temporary format)\n");
printf("%s %s %s %d %d %.4f %10.2f %.4f %s\n",matrixname,opname,sw,times,nnz,tottime,mxxops,bpnz,msstr);
end
if nargin <= 0
# DGEMV benchmark
for o=1024:1024
#for o=1024:256:2048*2
m=rand(o);
v=linspace(1,1,o)';
tic();
for i=1:n; m*v; end
t=toc();
Mflops=n*2.0*nnz(m)/(10^6 * t);
dgemvmflops=Mflops;
printf("%d GEMV for order %d in %g secs, so %10f Mflops\n",n,o,t,n*2.0*o*o/(10^6 * t));
end
quit
endif
# if nargin > 0, we continue
want_sparsersb_io=1;
if want_sparsersb_io != 1
source("ext/mminfo.m");
source("ext/mmread.m");
source("ext/mmwrite.m");
end
#matrices=ls("*.mtx")';
f=1;
uc=2;
while f<=nargin
MB=1024*1024;
mmn=cell2mat(argv()(f))';
mn=strtrim(mmn');
tic();
#nm=mmread(mn);
if want_sparsersb_io == 1
[nm,nrows,ncols,entries,rep,field,symm]=sparsersb(mn);
nm=sparse(nm);
if (symm=='S')
uc+=1;
end
else
[nm,nrows,ncols,entries,rep,field,symm]=mmread(mn);
#if(symm=="symmetric")uc+=2;endif
if(strcmp(symm,"symmetric"))uc+=1;endif
end
wr=0 ;
if wr==1
sparsersb(sparsersb(nm),"render",[mn,"-original.eps"]);
pct=-time;
#p=colamd(nm);
p=colperm(nm);
pct+=time;
pat=-time;
nm=nm(:,p);
pat+=time;
#sparsersb(sparsersb(nm),"render",[mn,"-colamd.eps"])
sparsersb(sparsersb(nm),"render",[mn,"-colperm.eps"]);
end
fsz=stat(mn).size;
rt=toc();
[ia,ja,va]=find(nm);
printf("%s: %.2f MBytes read in %.4f s (%10.2f MB/s)\n",mn',fsz/MB,rt,fsz/(rt*MB));
#ia=ia'; ja=ja'; va=va';
sep=" ";
csvlstring=sprintf("#mn entries nrows ncols");
csvdstring=sprintf("%%:%s%s%d%s%d%s%d",mn,sep,entries,sep,nrows,sep,ncols);
for ski=1:uc
oppnz=1;
# FIXME: what about symmetry ?
sparsekw="sparse";
if(ski==2)sparsekw="sparsersb";endif
if(ski==3);
oppnz=2;
sparsekw="sparsersb";
tic(); [nm]=sparsersb(mn); rt=toc();
sparsersb(nm,"info")
printf("%s: %.2f MBytes read by librsb in %.4f s (%10.2f MB/s)\n",mn',fsz/MB,rt,fsz/(rt*MB));
endif
if(ski==4);
nm=tril(nm);
endif
[ia,ja,va]=find(nm);
rnz=nnz(nm);
printf("benchmarking %s\n",sparsekw);
#printf("symmetry ? %s\n",issymmetric(sparse(nm)));
mrc=rows(nm); mcc=columns(nm);
if(ski!=3);
tic();
eval(["for i=1:n; om=",sparsekw,"(ia,ja,va,mrc,mcc,\"summation\"); end"]);
printf("benchmarking %s\n",sparsekw);
at=toc();
#if(ski==2) tic(); nm=sparsersb(om,"autotune","N");om=nm; att=toc(); ;endif
mnz=nnz(om);
amflops=n*2.0*mnz/(10^6 * at);
printf("%s (%s) %d spBLD for %d nnz in %.4f secs, so %10.2f Mflops\n",mn',sparsekw,n,rnz,at,amflops);
else
mnz=rnz;
end
#rm=sparsersb(ia,ja,va);# UNFINISHED
r=linspace(1,1,size(om,1))';
v=linspace(1,1,size(om,2))';
tic(); for i=1:n; r+=om *v; end; umt=toc();
UMflops=oppnz*n*2.0*mnz/(10^6 * umt);
printf("%s (%s) %d spMV for %d nnz in %.4f secs, so %10.2f Mflops\n",mn',sparsekw,n,mnz,umt, UMflops);
bpnz=-1; # FIXME: bytes per nonzero!
msstr="?";# FIXME: matrix structure string!
# FIXME: finish the following!
#printbenchline(mn',"spMV",sparsekw,n,mnz,umt, UMflops,bpnz,msstr);
#
tmp=r;r=v;v=tmp;
tic(); for i=1:n; r+=om.'*v; end; tmt=toc();
TMflops=oppnz*n*2.0*mnz/(10^6 * tmt);
printf("%s (%s) %d spMVT for %d nnz in %.4f secs, so %10.2f Mflops\n",mn',sparsekw,n,mnz,tmt, TMflops);
if(ski<3);
csvlstring=sprintf("%s%s",csvlstring," n at amflops umt UMflops tmt TMflops");
csvdstring=sprintf("%s%s%d%s%f%s%f%s%f%s%f%s%f%s%f",csvdstring,sep,n,sep,at,sep,amflops,sep,umt,sep,UMflops,sep,tmt,sep,TMflops);
endif
end
++f;
printf("%s\n",csvlstring);
printf("%s\n",csvdstring);
end
printf("benchmark terminated\n");
%%MatrixMarket matrix coordinate real general
% a positive definitive matrix, as in
% http://www.ncsa.uiuc.edu/UserInfo/Resources/Hardware/IBMp690/IBM/usr/lpp/essl.html.en_US/html/essl43.html
% * *
% | 99 12 13 14 15 16 |
% | 12 99 12 13 14 15 |
% | 13 12 99 12 13 14 |
% | 14 13 12 99 12 13 |
% | 15 14 13 12 99 12 |
% | 16 15 14 13 12 99 |
% * *
6 6 36
1 1 99
1 2 12
1 3 13
1 4 14
1 5 15
1 6 16
2 1 12
2 2 99
2 3 12
2 4 13
2 5 14
2 6 15
3 1 13
3 2 12
3 3 99
3 4 12
3 5 13
3 6 14
4 1 14
4 2 13
4 3 12
4 4 99
4 5 12
4 6 13
5 1 15
5 2 14
5 3 13
5 4 12
5 5 99
5 6 12
6 1 16
6 2 15
6 3 14
6 4 13
6 5 12
6 6 99
#
# Copyright (C) 2011-2015 Michele Martone <michelemartone _AT_ users.sourceforge.net>
#
# This program is free software; you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation; either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program; if not, see <http://www.gnu.org/licenses/>.
#
#
# This program shall attempt solution of a problem saved in the MATLAB format as for the University of Florida collection.
#
# e.g.: http://www.cise.ufl.edu/research/sparse/mat/Hamm/memplus.mat
#
# s=load("~/memplus.mat");
1; # This is a script
s=load(argv(){length(argv())});
n=rows(s.Problem.A);
minres=1e-7;
maxit = n;
#maxit = 100;
b=s.Problem.b;
#A=sparse(s.Problem.A);
A=sparsersb(s.Problem.A);
X0=[];
RELRES=2*minres;
TOTITER=0;
M1=[]; M2=[];
M1=spdiag(A)\ones(n,1);
M2=spdiag(ones(n,1));
while RELRES >= minres ;
tic; [X1, FLAG, RELRES, ITER] = pcg (A, b, minres, maxit, M1,M2,X0); odt=toc;
RELRES
ITER;
TOTITER=TOTITER+ITER
toc
X0=X1;
end
================================================================================
This (sparsersb) is a plugin to interface the "librsb" high performance sparse
matrix library to GNU Octave.
Author: Michele MARTONE
================================================================================
Build / use instructions (using pkg):
To use the configure auto-detected librsb from within Octave:
> pkg -local -verbose install sparsersb-1.0.1.tar.gz
> pkg load sparsersb
> help sparsersb
Alternatively:
tar czf sparsersb-1.0.1.tar.gz
cd sparsersb-1.0.1/src
./configure
make
make check
It is possible to provide the configure script with a librsb sources
archive you have downloaded separately and make it build on the fly for you
and use it; e.g.:
./configure --with-librsb-tarball=$HOME/librsb-1.2.0.tar.gz
or setting:
export LIBRSB_TARBALL=$HOME/librsb-1.2.0.tar.gz
before entering in Octave and building with pkg:
> pkg -local -verbose install sparsersb-1.0.1.tar.gz
On many systems, you will have to build librsb with the PIC (-fPIC on GCC)
option or you will get link time problems.
More configure options:
./configure --help
Usage instructions without using pkg:
# go to the directory where sparsersb.oct is located and run Octave:
octave
# you can use the sparsersb function, starting with e.g.:
> help sparsersb
Check out http://librsb.sf.net for the latest librsb.
================================================================================
-- Loadable Function: S = sparsersb (A)
-- Loadable Function: S = sparsersb (I, J, SV, M, N)
-- Loadable Function: S = sparsersb (I, J, SV, M, N, NZMAX)
-- Loadable Function: S = sparsersb (I, J, SV)
-- Loadable Function: S = sparsersb (M, N)
-- Loadable Function: S = sparsersb (I, J, SV, M, N, "unique")
-- Loadable Function: sparsersb ("set", OPN, OPV)
-- Loadable Function: V = sparsersb (S, "get", MIF)
-- Loadable Function: V = sparsersb (S, QS)
-- Loadable Function: sparsersb (A,"save",MTXFILENAME)
-- Loadable Function: [S, NROWS, NCOLS, NNZ, REPINFO, FIELD, SYMMETRY]
= sparsersb (MTXFILENAME, MTXTYPESTRING)
-- Loadable Function: sparsersb (S,"render", FILENAME[, RWIDTH,
RHEIGHT])
-- Loadable Function: [O =] sparsersb (S,"autotune"[, TRANSA, NRHS,
MAXR, TMAX, TN, SF])
Create or manipulate sparse matrices using the RSB format provided
by librsb, as similarly as possible to `sparse'.
If A is a full matrix, convert it to a sparse matrix
representation, removing all zero values in the process.
Given the integer index vectors I and J, and a 1-by-`nnz' vector
of real or complex values SV, construct the sparse matrix
`S(I(K),J(K)) = SV(K)' with overall dimensions M and N.
The argument `NZMAX' is ignored but accepted for compatibility
with MATLAB and `sparsersb'.
If M or N are not specified their values are derived from the
maximum index in the vectors I and J as given by `M = max (I)', `N
= max (J)'.
Can load a matrix from a Matrix Market matrix file named
MTXFILENAME. The optional argument MTXTYPESTRING can specify
either real ("D") or complex ("Z") type. Default is real. In the
case MTXFILENAME is "?", a string listing the available numerical
types with BLAS-style characters will be returned. If the file
turns out to contain a Matrix Market dense vector, this will be
loaded.
If "save" is specified, saves a sparse RSB matrix as a Matrix
Market matrix file named MTXFILENAME.
*Note*: if multiple values are specified with the same I, J
indices, the corresponding values in SV will be added.
The following are all equivalent:
s = sparsersb (i, j, s, m, n)
s = sparsersb (i, j, s, m, n, "summation")
s = sparsersb (i, j, s, m, n, "sum")
If the optional "unique" keyword is specified, then if more than
two values are specified for the same I, J indices, only the last
value will be used.
`sparsersb (M, N)' will create an empty MxN sparse matrix and is
equivalent to `sparsersb ([], [], [], M, N)'.
If M or N are not specified, then `M = max (I)', `N = max (J)'.
If OPN is a string representing a valid librsb option name and OPV
is a string representing a valid librsb option value, these will
be passed to the `rsb_lib_set_opt_str()' function.
If MIF is a string specifying a valid librsb matrix info string
(valid for librsb's `rsb_mtx_get_info_from_string()'), then the
corresponding value will be returned for matrix `A', in string
`V'. If MIF is the an empty string (""), matrix structure
information will be returned.
If A is a sparsersb matrix and QS is a string, QS will be
interpreted as a query string about matrix A. String `V' will be
returned. See librsb's `rsb_mtx_get_info_str()'. *Note*: this
feature is still incomplete, and whatever the value of QS, a
general information string will be returned.
If S is a sparsersb matrix and the "render" keyword is specified,
and FILENAME is a string, A will be rendered as an Encapsulated
Postscript file FILENAME. Optionally, width and height can be
specified in `RWIDTH, RHEIGHT'. Defaults are 512.
If S is a sparsersb matrix and the "autotune" keyword is
specified, autotuning of the matrix will take place, with SpMV and
autotuning parameters. After the "autotune" string, the remaining
parameters are optional. Parameter TRANSA specifies whether to
tune for untransposed ("N") or transposed ("T"); NRHS the number
of right hand sides; MAXR the number of tuning rounds; TMAX the
threads to use. If giving an output argument O, that will be
assigned to the autotuned matrix, and the input one A will remain
unchanged. See librsb documentation for `rsb_tune_spmm' to learn
more.
Please note that on `sparsersb' type variables are available most,
but not all of the operators available for `full' or `sparse'
typed variables.
See also: full, sparse
Additional help for built-in functions and operators is
available in the on-line version of the manual. Use the command
`doc <topic>' to search the manual index.
Help and information about Octave is also available on the WWW
at http://www.octave.org and via the help@octave.org
mailing list.
#
# Copyright (C) 2011-2015 Michele Martone <michelemartone _AT_ users.sourceforge.net>
#
# This program is free software; you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation; either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program; if not, see <http://www.gnu.org/licenses/>.
#
# Undocumented
# function mbench(fname)
function matlabbench(fname)
addpath ext/
mp=matlabpath();
n=10;
for f=1:nargin
MB=1024*1024;
mmn=fname';
mn=strtrim(mmn');
tic();
matlabpath('ext');
mn;
nm=mmreadm(mn);
matlabpath(mp);
%fsz=stat(mn).size;
rnz=nnz(nm);
rt=toc();
[ia,ja,va]=find(nm);
%printf('%s: %.2f MBytes read in %.4f s (%10.2f MB/s)\n',mn',fsz/MB,rt,fsz/(rt*MB));
%ia=ia'; ja=ja'; va=va';
for ski=1:1
% FIXME: what about symmetry ?
sparsekw='sparse';
if(ski==2)sparsekw='sparsersb';end
tic();
for i=1:n; om=sparse(ia,ja,va); end
at=toc();
mnz=nnz(om);
amflops=n*2.0*mnz/(10^6 * at);
%printf('%s (%s) %d spBLD for %d nnz in %.4f secs, so %10.2f Mflops\n',mn',sparsekw,n,rnz,at,amflops);
amflops
%rm=sparsersb(ia,ja,va);% UNFINISHED
v=linspace(1,1,size(om,1))';
r=v;
tic(); for i=1:n r=r+om*v; end ; umt=toc();
UMflops=n*2.0*mnz/(10^6 * umt);
UMflops
%printf('%s (%s) %d spMV for %d nnz in %.4f secs, so %10.2f Mflops\n',mn',sparsekw,n,mnz,umt, UMflops);
tic(); for i=1:n r=r+om.'*v; end ; tmt=toc();
TMflops=n*2.0*mnz/(10^6 * tmt);
TMflops
%printf('%s (%s) %d spMV for %d nnz in %.4f secs, so %10.2f Mflops\n',mn',sparsekw,n,mnz,tmt, TMflops);
end
end
%printf('benchmark terminated successfully\n');
quit
end
#!/usr/bin/octave -q
#
# Copyright (C) 2011-2016 Michele Martone <michelemartone _AT_ users.sourceforge.net>
#
# This program is free software; you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation; either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program; if not, see <http://www.gnu.org/licenses/>.
#
# TODO: document this file's functions so they get installed and are properly usable.
# TODO: sprand should not be used in a consistent way
1; # This is a script.
pkg load sparsersb
function dt=sparsersbbench__(precmd,cmd,postcmd,mint)
# ..
eval(precmd);
nops=0;
tic();
do
++nops;
eval(cmd);