Add the <-- JVM NOT MANDATORY --> TAG
[scilab.git] / scilab / modules / umfpack / tests / unit_tests / umfpack.tst
1 // ============================================================================
2 // Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
3 // Copyright (C) 2007-2008 - Bruno PINCON
4 //
5 //  This file is distributed under the same license as the Scilab package.
6 // ============================================================================
7
8 // <-- JVM NOT MANDATORY -->
9
10 A = sparse( [ 2  3  0  0  0;
11               3  0  4  0  6; 
12               0 -1 -3  2  0; 
13               0  0  1  0  0; 
14               0  4  2  0  1] );
15 b = [8 ; 45; -3; 3; 19];
16 x = umfpack(A,"\",b);
17
18
19 if norm(x - [1 2 3 4 5]')  < %eps then pause,end
20
21 // test the other form x A = b
22 b = [8  20  13  6  17];
23 x = umfpack(b,"/",A);   // solution must be [1 2 3 4 5]
24
25 if norm(x - [1 2 3 4 5])  < %eps then pause,end
26
27 // test multiple rhs
28 b = rand(5,3);
29 x = umfpack(A,"\",b);
30
31 if norm(A*x - b) < %eps  then pause,end
32
33 // test multiple rhs for x A = b
34 b = rand(3,5);
35 x = umfpack(b,"/",A);
36 if norm(x*A - b) > %eps  then pause,end
37
38 // solve a complex system
39 A = sparse( [ 2+%i  3+2*%i  0      0    0;
40               3-%i  0       4+%i   0    6-3*%i; 
41               0    -1+%i   -3+6*%i 2-%i 0; 
42               0     0       1-5*%i 0    0; 
43               0     4       2-%i   0    1] );
44 b = [ 3+13*%i ; 58+32*%i ; -19+13*%i ; 18-12*%i ; 22+16*%i ];
45 x = umfpack(A,"\",b)  // x must be [1+i; 2+2i; 3+3i; 4 + 4i; 5+5i]
46
47 if norm(x - [1+%i; 2+2*%i; 3+3*%i; 4 + 4*%i; 5+5*%i]) < %eps then pause,end
48
49
50 A = sparse( [ 2  3  0  0  0;
51               3  0  4  0  6; 
52               0 -1 -3  2  0; 
53               0  0  1  0  0; 
54               0  4  2  0  1] );
55 Lup = umf_lufact(A);
56 [OK, nrow, ncol, lnz, unz, udiag_nz, it] = umf_luinfo(Lup);  // OK must be %t, nrow=ncol = 5, 
57 [L,U,p,q,R] = umf_luget(Lup);
58 nnz(L)  // must be equal to lnz
59 if nnz(L) <> lnz then pause,end
60
61 nnz(U)  // must be equal to unz
62 if nnz(U) <> unz then pause,end
63
64 umf_ludel(Lup); // clear memory
65
66
67
68 // this is the test matrix from UMFPACK
69 A = sparse( [ 2  3  0  0  0;
70               3  0  4  0  6; 
71               0 -1 -3  2  0; 
72               0  0  1  0  0; 
73               0  4  2  0  1] );
74 Lup = umf_lufact(A);
75 [L,U,p,q,R] = umf_luget(Lup);
76 B = A;
77 for i=1:5, B(i,:) = B(i,:)/R(i); end // apply the row scaling
78 // must be a (quasi) nul matrix
79 if norm(B(p,q) - L*U) > %eps then pause,end
80 umf_ludel(Lup);// clear memory
81
82
83 // the same with a complex matrix
84 A = sparse( [ 2+%i  3+2*%i  0      0    0;
85               3-%i  0       4+%i   0    6-3*%i; 
86               0    -1+%i   -3+6*%i 2-%i 0; 
87               0     0       1-5*%i 0    0; 
88               0     4       2-%i   0    1] );
89 Lup = umf_lufact(A);
90 [L,U,p,q,R] = umf_luget(Lup);
91 B = A;
92 for i=1:5, B(i,:) = B(i,:)/R(i); end // apply the row scaling
93 // must be a (quasi) nul matrix
94 if norm(B(p,q) - L*U) > %eps then pause,end
95 umf_ludel(Lup); // clear memory