-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathdriv_test.m
66 lines (55 loc) · 2.08 KB
/
driv_test.m
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
numsweeps = 7; % useful cases: 1 and small values and large values
numthreads = 4; % useful cases: 1 and 4 and larger values
a = delsq(numgrid('S',40));
a = diagscale(a);
u_ichol = ichol(a)';
[l_ilu u_ilu] = ilu(a);
% test paric_ref
u = triu(a);
u = paric_ref(a, u, numsweeps);
fprintf('paric_ref: %e\n', norm(u-u_ichol,'fro')/norm(u_ichol,'fro'));
% test paric_ref with residual history computation
u = triu(a);
[u resid] = paric_ref(a, u, numsweeps);
fprintf('paric_ref: %e\n', norm(u-u_ichol,'fro')/norm(u_ichol,'fro'));
% test paric/synchronous
u = triu(a);
u = paric(a, u, 1, numsweeps, numthreads);
fprintf('paric/sync: %e\n', norm(u-u_ichol,'fro')/norm(u_ichol,'fro'));
% test paric/asynchronous
u = triu(a);
u = paric(a, u, 0, numsweeps, numthreads);
fprintf('paric/async: %e\n', norm(u-u_ichol,'fro')/norm(u_ichol,'fro'));
% test parilu_ref
l = tril(a);
u = triu(a);
[l u] = parilu_ref(a, l, u, numsweeps);
fprintf('parilu_ref:L %e\n', norm(l-l_ilu,'fro')/norm(l_ilu,'fro'));
fprintf('parilu_ref:U %e\n', norm(u-u_ilu,'fro')/norm(u_ilu,'fro'));
% test parilu/synchronous
l = tril(a);
u = triu(a);
[l u] = parilu(a, l, u, 1, numsweeps, numthreads);
fprintf('parilu/sync:L %e\n', norm(l-l_ilu,'fro')/norm(l_ilu,'fro'));
fprintf('parilu/sync:U %e\n', norm(u-u_ilu,'fro')/norm(u_ilu,'fro'));
% test parilu/asynchronous
l = tril(a);
u = triu(a);
[l u] = parilu(a, l, u, 0, numsweeps, numthreads);
fprintf('parilu/async:L %e\n', norm(l-l_ilu,'fro')/norm(l_ilu,'fro'));
fprintf('parilu/async:U %e\n', norm(u-u_ilu,'fro')/norm(u_ilu,'fro'));
% test pariutdt_ref (sync)
u = triu(a);
d = speye(length(a));
[u d resid] = pariutdu_ref(a, u, d, numsweeps);
fprintf('pariutdu_ref: %e\n', norm(u-l_ilu','fro')/norm(l_ilu,'fro'));
% test pariutdt/synchronous
u = triu(a);
d = speye(length(a));
[u d] = pariutdu(a, u, d, 1, numsweeps, numthreads);
fprintf('pariutdu/sync: %e\n', norm(u-l_ilu','fro')/norm(l_ilu,'fro'));
% test pariutdt/asynchronous
u = triu(a);
d = speye(length(a));
[u d] = pariutdu(a, u, d, 0, numsweeps, numthreads);
fprintf('pariutdu/asyn: %e\n', norm(u-l_ilu','fro')/norm(l_ilu,'fro'));