1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
|
% ********************************************************************
% Copyright (C) 2023 Advanced Micro Devices, Inc. All rights reserved.
% ********************************************************************
function genmat_posdef(n,nz)
%{
-----------------------------------------------------------
This function generates all the files required for a
test case of the Cholesky re-factorization functionality in rocSOLVER.
genmat_gen(n, nz)
Inputs:
n Size of the problem (number of row and columns
of the sparse positive-definite matrix A)
nz Max number of non-zero elements in A
Outputs:
isok This will be -1 if any of the files could not
be generated
The resulting files are created in the working directory.
See printall.m for more details.
(This script is for internal use only. It is not part of
rocSOLVER library interface and could change or be removed
without any notice)
-----------------------------------------------------------
%}
%---------------------------
% generate random sparse matrix
% with nz non-zero elements
%---------------------------
per = nz / (n*n);
A = sprand(n, n, per);
%---------------------------
% Take lower part L and make it diagonally dominant
% (this means that the symmetric matrix L+L^T will be
% positive definite)
%---------------------------
A=tril(A);
A = fixmat(A, n/2);
%---------------------------
% If necessary, remove extra non-zero elements
% so that the final symmetric matrix has no more
% than nz non-zero elements
%---------------------------
nzsym = (nz + n) / 2;
dif = nnz(A) - nzsym;
if(dif > 0),
A = cleanmat(A, dif);
end;
%---------------------------
% compute the final symmetric matrix
% (functions work with only the lower part, but better
% to have a full symmetric to generate the linear systems)
%---------------------------
A = A + transpose(A);
%---------------------------
% Compute Cholesky factorization
%---------------------------
[T, p, Q] = chol(A, 'vector', 'lower');
%---------------------------
% Generate output files
%---------------------------
disp('')
disp('Generating files...')
isok = printall("po", A, T, Q);
if (isok),
disp('OK')
disp('')
else
disp('FAILURE!');
disp('')
end;
%------------------------
% quick validation
%------------------------
ptr = dlmread('ptrT');
ind = dlmread('indT');
val = dlmread('valT');
t = csr2full(n, ptr, ind, val);
qq = dlmread('Q') + 1;
q = eye(n,n);
q = q(:,qq);
disp('Quick checking:')
b = dlmread('B_1');
x = dlmread('X_1');
err1 = norm(q*t*t'*q'*x-b, "fro")
b = dlmread('B_10');
x = dlmread('X_10');
err2 = norm(q*t*t'*q'*x-b, "fro")
b = dlmread('B_30');
x = dlmread('X_30');
err3 = norm(q*t*t'*q'*x-b, "fro")
disp('')
end
|
