Added the .dia.ref files back in the unit tests.
[scilab.git] / scilab / modules / linear_algebra / tests / unit_tests / bdiag.dia.ref
1 // =============================================================================
2 // Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
3 // Copyright (C) ????-2008 - INRIA Michael Baudin
4 //
5 //  This file is distributed under the same license as the Scilab package.
6 // =============================================================================
7 function r=Err(x),r=norm(x,1),endfunction
8 rand('normal')
9 //==========================================================================
10 //==============================   bdiag      ==============================
11 //==========================================================================
12 if bdiag([])<>[] then bugmes();quit;end
13 [ab,x]=bdiag([]);
14 if ab<>[]|x<>[] then bugmes();quit;end
15 [ab,x,bs]=bdiag([]);
16 if ab<>[]|x<>[]|bs<>[] then bugmes();quit;end
17 if execstr('bdiag([1 2;3 4;5 6])','errcatch')==0 then bugmes();quit;end
18 //Small dimension
19 //---------------
20 //Real case
21 e=1.d-1;
22 A=[1 1  2 3 4 5
23    0 1  6 7 8 9
24    0 0  1 e 3 1
25    0 0 -e 1 5 9
26    0 0  0 0 2 e
27    0 0  0 0 0 3];
28 X1=[0.5,0.3,0,0.3,0.3,0.2;
29    1,0.6,0.5,0.1,0.7,0.4;
30    0.7,0.1,0.4,0.6,0.1,1;
31    0,0.6,0.2,0.3,0.4,0.5;
32    0.6,0.7,0.5,0.7,0.7,0.5;
33    0.3,0.3,0.4,0.5,0.9,0.6]
34  X1  =
35  
36     0.5    0.3    0.     0.3    0.3    0.2  
37     1.     0.6    0.5    0.1    0.7    0.4  
38     0.7    0.1    0.4    0.6    0.1    1.   
39     0.     0.6    0.2    0.3    0.4    0.5  
40     0.6    0.7    0.5    0.7    0.7    0.5  
41     0.3    0.3    0.4    0.5    0.9    0.6  
42 A=inv(X1)*A*X1;
43 Ab1=bdiag(A);
44 if or(triu(Ab1,-1)<>Ab1) then bugmes();quit;end
45 [Ab2,X]=bdiag(A);
46 if Err(Ab2-Ab1)>>10*%eps then bugmes();quit;end
47 if Err(Ab2-inv(X)*A*X )>1d6*%eps then bugmes();quit;end
48 [Ab2,X,bs]=bdiag(A);
49 if Err(Ab2-Ab1)>>10*%eps then bugmes();quit;end
50 if Err(Ab2-inv(X)*A*X )>2.d-10 then bugmes();quit;end
51 if or(size(bs)<>[3,1]) then bugmes();quit;end
52 if sum(bs)<>size(A,1) then bugmes();quit;end
53 if or(bs<=0) then bugmes();quit;end
54 [Ab2,X,bs]=bdiag(A,1);
55 if Err(Ab2-Ab1)>>10*%eps then bugmes();quit;end
56 if Err(Ab2-inv(X)*A*X )>2d-7 then bugmes();quit;end
57 if or(size(bs)<>[1,1]) then bugmes();quit;end
58 if sum(bs)<>size(A,1) then bugmes();quit;end
59 if or(bs<=0) then bugmes();quit;end
60 //Complex case
61 e=1.d-1;
62 A=[1 1  2 3 4 5
63    0 1  6 7 8 9
64    0 0  1 e 3 1
65    0 0 -e 1 5 9
66    0 0  0 0 2 e
67    0 0  0 0 0 3];
68 X1=[0.5,0.3,0,0.3,0.3,0.2;
69    1,0.6,0.5,0.1,0.7,0.4;
70    0.7,0.1,0.4,0.6,0.1,1;
71    0,0.6,0.2,0.3,0.4,0.5;
72    0.6,0.7,0.5,0.7,0.7,0.5;
73    0.3,0.3,0.4,0.5,0.9,0.6]+%i*eye(A);
74 A=inv(X1)*A*X1;
75 Ab1=bdiag(A);
76 if or(triu(Ab1)<>Ab1) then bugmes();quit;end
77 [Ab2,X]=bdiag(A);
78 if Err(Ab2-Ab1)>>10*%eps then bugmes();quit;end
79 if Err(Ab2-inv(X)*A*X )>1.d-8 then bugmes();quit;end
80 [Ab2,X,bs]=bdiag(A);
81 if Err(Ab2-Ab1)>>10*%eps then bugmes();quit;end
82 if Err(Ab2-inv(X)*A*X )>1.d-8 then bugmes();quit;end
83 if size(bs,2)<>1 then bugmes();quit;end
84 if sum(bs)<>size(A,1) then bugmes();quit;end
85 if or(bs<=0) then bugmes();quit;end
86 //Large dimension
87 //---------------
88 //Real case
89 A=rand(25,25);
90 Ab1=bdiag(A);
91 if or(triu(Ab1,-1)<>Ab1) then bugmes();quit;end
92 [Ab2,X]=bdiag(A);
93 if Err(Ab2-Ab1)>>10*%eps then bugmes();quit;end
94 if Err(Ab2-inv(X)*A*X )>10000*%eps then bugmes();quit;end
95 [Ab2,X,bs]=bdiag(A);
96 if Err(Ab2-Ab1)>>10*%eps then bugmes();quit;end
97 if Err(Ab2-inv(X)*A*X )>10000*%eps then bugmes();quit;end
98 if size(bs,2)<>1 then bugmes();quit;end
99 if sum(bs)<>size(A,1) then bugmes();quit;end
100 if or(bs<=0) then bugmes();quit;end
101 //Complex case
102 A=rand(25,25)+%i*rand(25,25);
103 Ab1=bdiag(A);
104 if or(triu(Ab1)<>Ab1) then bugmes();quit;end
105 [Ab2,X]=bdiag(A);
106 if Err(Ab2-Ab1)>>10*%eps then bugmes();quit;end
107 if Err(Ab2-inv(X)*A*X )>10000*%eps then bugmes();quit;end
108 [Ab2,X,bs]=bdiag(A);
109 if Err(Ab2-Ab1)>>10*%eps then bugmes();quit;end
110 if Err(Ab2-inv(X)*A*X )>10000*%eps then bugmes();quit;end
111 if size(bs,2)<>1 then bugmes();quit;end
112 if sum(bs)<>size(A,1) then bugmes();quit;end
113 if or(bs<=0) then bugmes();quit;end