1function ss_index = ssindex (matrixlist)
2%SSINDEX create the index for the SuiteSparse Matrix Collection
3%
4% ss_index = ssindex (matrixlist)
5%
6% matrixlist: an integer list, in the range of 1 to the length of the
7%   ss_listing.txt file, containing a list of matrices for which to modify
8%   the ss_index entries.  If matrixlist is not present, then the ss_index
9%   is created from scratch.
10%
11% ss_index:  a struct containing the index information, with the
12%   following fields, assuming that there are n matrices in the collection:
13%
14%   LastRevisionDate: a string with the date and time the index was updated.
15%   DowloadTimeStamp: date and time the index was last downloaded.
16%
17%   Group: an n-by-1 cell array.  Group {i} is the group for matrix i.
18%   Name: n-by-1 cell array.  Name {i} is the name of matrix i.
19%
20%   The following fields are n-by-1 vectors unless otherwise specified.
21%   nrows(id) gives the number of rows of the matrix with id = Problem.id,
22%   for example.
23%
24%   nrows           number of rows
25%   ncols           number of columns
26%   nnz             number of entries in A
27%   RBtype          Rutherford/Boeing type, an n-by-3 char array
28%   isBinary        1 if binary, 0 otherwise
29%   isReal          1 if real, 0 if complex
30%   cholcand        1 if a candidate for sparse Cholesky, 0 otherwise
31%   numerical_symmetry  numeric symmetry (0 to 1, where 1=symmetric)
32%   pattern_symmetry    pattern symmetry (0 to 1, where 1=symmetric)
33%   nnzdiag         nnz (diag (A)) if A is square, 0 otherwise
34%   nzero           nnz (Problem.Zeros)
35%   nentries        nnz + nzero
36%   amd_lnz         nnz(L) for chol(C(p,p)) where, C=A+A', p=amd(C)
37%   amd_flops       flop count for chol(C(p,p)) where, C=A+A', p=amd(C)
38%   amd_vnz         nnz in Householder vectors for qr(A(:,colamd(A)))
39%   amd_rnz         nnz in R for qr(A(:,colamd(A)))
40%   nblocks         # of blocks from dmperm
41%   sprank          sprank(A)
42%   ncc             # of strongly connected components
43%   posdef          1 if positive definite, 0 otherwise
44%   isND            1 if a 2D/3D problem, 0 otherwise
45%   isGraph         1 if a graph, 0 otherwise
46%   lowerbandwidth      lower bandwidth, [i j]=find(A), max(0,max(i-j))
47%   upperbandwidth      upper bandwidth, [i j]=find(A), max(0,max(j-i))
48%   rcm_lowerbandwidth  lower bandwidth after symrcm
49%   rcm_upperbandwidth  upper bandwidth after symrcm
50%   xmin            smallest nonzero value
51%   xmax            largest nonzero value
52%
53% If the statistic is intentionally not computed, it is set to -2.  Some
54% statistics are not computed for rectangular or structurally singular
55% matrices, for example.  If an attempt to compute the statistic was made, but
56% failed, it is set to -1.  If no attempt yet has been made to compute the
57% entry, it is set to -3.
58%
59% Example:
60%   ssindex
61%   ssindex (267:300)
62%
63% If updating the ss_index.mat file, the function first loads it from its
64% default location, via ssget.  This function then saves the new ss_index into
65% the ss_index.mat file in the current working directory (not overwriting the
66% old copy, unless it resides in the current working directory).  It creates
67% the ssstats.csv file used by ssgui.java and sskinds.m and places it in the
68% current working directory.
69%
70% See also ssstats.
71
72% Copyright 2006-2019, Timothy A. Davis
73
74% Requires the SuiteSparse set of packages: CHOLMOD, RBio, CSparse
75
76%   10/13/2001: Created by Erich Mirabal
77
78%-------------------------------------------------------------------------------
79% initialize an empty index
80%-------------------------------------------------------------------------------
81
82% load the filenames
83topdir = sslocation ;
84files = textread ([topdir 'files/ss_listing.txt'], '%s') ;
85
86% if no input, assume we have to do the whole file list
87create_new = 0 ;
88if (nargin < 1)
89    matrixlist = 1:length(files) ;
90    create_new = 1 ;
91elseif (~isempty (matrixlist))
92    % validate the input : range is limited by the files variable
93    if (min (matrixlist) < 1) || (max (matrixlist) > length (files))
94        error ('%s: %s', mfilename, 'Invalid input parameter.') ;
95    end
96end
97
98if (~create_new)
99    % load the index from file
100    fprintf ('Loading existing ss_index.mat file\n') ;
101    try
102        load ss_index
103        fprintf ('loaded ss_index in current directory:\n%s\n', pwd) ;
104        dir
105    catch
106        fprintf ('loading ss_index = ssget\n') ;
107        ss_index = ssget ;
108    end
109end
110
111% get the current sskinds
112if (create_new)
113    kinds = cell (length (files),1) ;
114else
115    try
116        kinds = sskinds ;
117    catch
118        kinds = cell (2,1) ;
119    end
120    for i = matrixlist
121        kinds {i} = '' ;
122    end
123end
124
125% revision tracking device
126ss_index.LastRevisionDate = datestr (now) ;
127
128% the index structure needs a download date for version tracking
129ss_index.DownloadTimeStamp = now ;
130
131% start the index from scratch
132if (create_new)
133
134    fprintf ('Creating new ss_index.mat file\n') ;
135    nothing = -3 * ones (1, length (files)) ;
136
137    ss_index.Group = cell (size (files)) ;
138    ss_index.Name = cell (size (files)) ;
139    ss_index.nrows = nothing ;
140    ss_index.ncols = nothing ;
141    ss_index.nnz = nothing ;
142    ss_index.nzero = nothing ;
143    ss_index.nentries = nothing ;
144    ss_index.pattern_symmetry = nothing ;
145    ss_index.numerical_symmetry = nothing ;
146    ss_index.isBinary = nothing ;
147    ss_index.isReal = nothing ;
148    ss_index.nnzdiag = nothing ;
149    ss_index.posdef = nothing ;
150    ss_index.amd_lnz   = nothing ;
151    ss_index.amd_flops = nothing ;
152    ss_index.amd_vnz   = nothing ;
153    ss_index.amd_rnz   = nothing ;
154    ss_index.nblocks   = nothing ;
155    ss_index.sprank    = nothing ;
156    ss_index.RBtype = char (' '*ones (length (files),3)) ;
157    ss_index.cholcand = nothing ;
158    ss_index.ncc = nothing ;
159    ss_index.isND = nothing ;
160    ss_index.isGraph = nothing ;
161    ss_index.lowerbandwidth = nothing ;
162    ss_index.upperbandwidth = nothing ;
163    ss_index.rcm_lowerbandwidth = nothing ;
164    ss_index.rcm_upperbandwidth = nothing ;
165    ss_index.xmin = nothing ;
166    ss_index.xmax = nothing ;
167
168else
169
170    % make sure we have the right length for the arrays
171    if length (ss_index.nrows) < max (matrixlist)
172
173        len = max (matrixlist) - length (ss_index.nrows) ;
174        nothing = -ones (1, len) ;
175
176        if (len > 0)
177            for i = matrixlist
178                ss_index.Group {i} = '' ;
179                ss_index.Name {i} = '' ;
180            end
181            ss_index.nrows      = [ss_index.nrows nothing] ;
182            ss_index.ncols      = [ss_index.ncols nothing] ;
183            ss_index.nnz        = [ss_index.nnz nothing] ;
184            ss_index.nzero      = [ss_index.nzero nothing] ;
185            ss_index.nentries   = [ss_index.nentries nothing] ;
186            ss_index.pattern_symmetry = [ss_index.pattern_symmetry nothing] ;
187            ss_index.numerical_symmetry = [ss_index.numerical_symmetry nothing];
188            ss_index.isBinary   = [ss_index.isBinary nothing] ;
189            ss_index.isReal     = [ss_index.isReal nothing] ;
190            ss_index.nnzdiag    = [ss_index.nnzdiag nothing] ;
191            ss_index.posdef     = [ss_index.posdef nothing] ;
192            ss_index.amd_lnz    = [ss_index.amd_lnz nothing] ;
193            ss_index.amd_flops  = [ss_index.amd_flops nothing] ;
194            ss_index.amd_vnz    = [ss_index.amd_vnz nothing] ;
195            ss_index.amd_rnz    = [ss_index.amd_rnz nothing] ;
196            ss_index.nblocks    = [ss_index.nblocks nothing] ;
197            ss_index.sprank     = [ss_index.sprank nothing] ;
198            ss_index.RBtype     = [ss_index.RBtype ; char (' '*ones (len,3))] ;
199            ss_index.cholcand   = [ss_index.cholcand nothing] ;
200            ss_index.ncc        = [ss_index.ncc nothing] ;
201            ss_index.isND       = [ss_index.isND nothing] ;
202            ss_index.isGraph    = [ss_index.isGraph nothing] ;
203            ss_index.lowerbandwidth     = [ss_index.lowerbandwidth nothing] ;
204            ss_index.upperbandwidth     = [ss_index.upperbandwidth nothing] ;
205            ss_index.rcm_lowerbandwidth = [ss_index.rcm_upperbandwidth nothing];
206            ss_index.rcm_upperbandwidth = [ss_index.rcm_upperbandwidth nothing];
207            ss_index.xmin = [ss_index.xmin nothing] ;
208            ss_index.xmax = [ss_index.xmax nothing] ;
209        end
210    end
211end
212
213if (length (matrixlist) > 0)
214    fprintf ('Will process %d files\n', length (matrixlist)) ;
215end
216
217nmat = length (ss_index.nrows) ;
218filesize = zeros (nmat,1) ;
219
220%-------------------------------------------------------------------------------
221% look through the directory listing
222%-------------------------------------------------------------------------------
223
224for i = matrixlist
225
226    % note that the matrix is not loaded in this for loop
227    ffile = deblank (files {i}) ;
228
229    % group is the first part of the string up to the character before
230    % the last file separator
231    gi = find (ffile == '/') ;
232    gi = gi (end) ;
233    groupN = char (ffile (1:gi-1)) ;
234
235    % name is the last section of the string after the last file separator
236    matrixN = char (ffile (gi+1:end)) ;
237
238    % get the directory info of the .mat file
239    fileInfo = dir ([topdir 'mat/' ffile '.mat']) ;
240
241    % set the file's data into the data arrays
242    ss_index.Name {i} = matrixN ;
243    ss_index.Group {i} = groupN ;
244
245    if (length (fileInfo) > 0)                                              %#ok
246        filesize (i) = fileInfo.bytes ;
247    else
248        filesize (i) = -1 ;
249    end
250
251end
252
253if (length (matrixlist) > 0)
254    fprintf ('\n======================================================\n') ;
255    fprintf ('Matrices will processed in the following order:\n') ;
256    for i = matrixlist
257        ffile = deblank (files {i}) ;
258        fprintf ('Matrix %d: %s filesize %d\n', i, ffile, filesize (i)) ;
259        if (filesize (i) == -1)
260            fprintf ('skip this file (not found)\n') ;
261            continue ;
262        end
263    end
264end
265
266fprintf ('Hit enter to continue\n') ;
267pause
268
269%-------------------------------------------------------------------------------
270% load the matrices
271%-------------------------------------------------------------------------------
272
273% known to be positive definite / indefinite:
274known_posdef = [ 939 1252 1267 1268 1423 1453 1455 2541:2547 ] ;
275known_indef = [ 1348:1368 1586 1411 1901:1905] ;
276
277% known to be irreducible, but dmperm takes too long:
278known_irreducible = [ 1902:1905 ] ;
279% known_irreducible = [ ] ;
280
281t = tic ;
282
283for k = 1:length (matrixlist)
284
285    %---------------------------------------------------------------------------
286    % get the matrix
287    %---------------------------------------------------------------------------
288
289    id = matrixlist (k) ;
290    ffile = deblank (files {id}) ;
291    fprintf ('\n============================== Matrix %d: %s\n', id, ffile) ;
292    if (filesize (id) == -1)
293	fprintf ('skip this file\n') ;
294	continue ;
295    end
296    load ([topdir 'mat/' ffile]) ;
297
298    % display the Problem struct
299    disp (Problem) ;
300
301    %---------------------------------------------------------------------------
302    % get all stats
303    %---------------------------------------------------------------------------
304
305    kinds {id} = Problem.kind ;
306
307    fprintf ('%s/%s\n', ss_index.Group {id}, ss_index.Name {id}) ;
308
309    if (~isequal (Problem.name, [ss_index.Group{id} '/' ss_index.Name{id}]))
310        error ('name mismatch!') ;
311    end
312    if (Problem.id ~= id)
313        error ('id mismatch!') ;
314    end
315
316    skip_chol = (any (id == known_posdef) || any (id == known_indef)) ;
317    skip_dmperm = any (id == known_irreducible) ;
318
319    if (isfield (Problem, 'Zeros'))
320	stats = ssstats (Problem.A, Problem.kind, skip_chol, ...
321            skip_dmperm, Problem.Zeros) ;
322    else
323	stats = ssstats (Problem.A, Problem.kind, skip_chol, ...
324            skip_dmperm) ;
325    end
326
327    %---------------------------------------------------------------------------
328    % fix special cases
329    %---------------------------------------------------------------------------
330
331    if (stats.posdef < 0)
332	if (any (id == known_posdef))
333	    fprintf ('known posdef\n') ;
334	    stats.posdef = 1 ;
335	elseif (any (id == known_indef))
336	    fprintf ('known indef\n') ;
337	    stats.posdef = 0 ;
338	end
339    end
340
341    if (any (id == known_irreducible) && stats.sprank < 0)
342	% full sprank, and not reducible to block triangular form,
343	% but dmperm takes too long
344	fprintf ('known irreducible\n') ;
345	stats.sprank = stats.nrows  ;
346	stats.nblocks = 1 ;
347    end
348
349    % display the stats
350    disp (stats) ;
351
352    %---------------------------------------------------------------------------
353    % save the stats in the index
354    %---------------------------------------------------------------------------
355
356    ss_index.nrows (id) = stats.nrows ;
357    ss_index.ncols (id) = stats.ncols ;
358    ss_index.nnz (id) = stats.nnz ;
359    ss_index.nzero (id) = stats.nzero ;
360    ss_index.nentries (id) = stats.nentries ;
361    ss_index.pattern_symmetry (id) = stats.pattern_symmetry ;
362    ss_index.numerical_symmetry (id) = stats.numerical_symmetry ;
363    ss_index.isBinary (id) = stats.isBinary ;
364    ss_index.isReal (id) = stats.isReal ;
365    ss_index.nnzdiag (id) = stats.nnzdiag ;
366    ss_index.posdef (id) = stats.posdef ;
367    ss_index.amd_lnz (id) = stats.amd_lnz ;
368    ss_index.amd_flops (id) = stats.amd_flops ;
369    ss_index.amd_vnz (id) = stats.amd_vnz ;
370    ss_index.amd_rnz (id) = stats.amd_rnz ;
371    ss_index.nblocks (id) = stats.nblocks ;
372    ss_index.sprank (id) = stats.sprank ;
373    ss_index.RBtype (id,:) = stats.RBtype ;
374    ss_index.cholcand (id) = stats.cholcand ;
375    ss_index.ncc (id) = stats.ncc ;
376    ss_index.isND (id) = stats.isND ;
377    ss_index.isGraph (id) = stats.isGraph ;
378    ss_index.lowerbandwidth (id) = stats.lowerbandwidth ;
379    ss_index.upperbandwidth (id) = stats.upperbandwidth ;
380    ss_index.rcm_lowerbandwidth (id) = stats.rcm_lowerbandwidth ;
381    ss_index.rcm_upperbandwidth (id) = stats.rcm_upperbandwidth ;
382    ss_index.xmin (id) = stats.xmin ;
383    ss_index.xmax (id) = stats.xmax ;
384
385    %---------------------------------------------------------------------------
386    % clear the problem and save the index and ssstats.csv
387    %---------------------------------------------------------------------------
388
389    clear Problem
390    fprintf ('time since last save: %g\n', toc (t)) ;
391    if (toc (t) > 20 || k == length (matrixlist))
392        t = tic ;
393        fprintf ('\n ... saving ss_index ...\n') ;
394        save ss_index ss_index
395
396        fprintf ('\nCreating ssstats.csv in current directory:\n')
397        fprintf ('%s/ssstats.csv\n', pwd) ;
398
399        sscsv_write (ss_index, kinds) ;
400
401        % flush the diary
402        if (strcmp (get (0, 'Diary'), 'on'))
403            diary off
404            diary on
405        end
406    end
407end
408
409