* Bug #8210 fixed - Inserting UMFPACK examples in the Demos gui as demos
[scilab.git] / scilab / modules / umfpack / demos / harwell-boeing.sce
1 // Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
2 // Copyright Bruno Pinçon, ESIAL-IECN, Inria CORIDA project
3 //   <bruno.pincon@iecn.u-nancy.fr>
4 //
5 // This file is released under the 3-clause BSD license. See COPYING-BSD.
6
7 function demos_harwell_boeing()
8
9     format(15);
10     cd "SCI/modules/umfpack/demos/"
11
12     disp(_("This test compares umfpack v3 and sparse v1.3 via their Scilab interface."));
13     disp(_("The test consists in loading Harwell-Boeing sparse matrices and solve linear system with a random rhs."));
14     disp(_("The matrices presented here come from the NIST server Matrix server:"));
15     printf("     http://math.nist.gov/MatrixMarket/matrices.html\n");
16     disp(_("Warning: tests 2 and 3 take much more time than the others."));
17     halt(_("\nPress Return to continue...\n"));
18
19
20     ////////////////////////////////////////////////////////////////////////////////////////////////////
21     // Test 1: Jacobian of a nonlinear system of ordinary differential equations (ODEs) modeling a laser
22     clear
23     nb_rep = 10;
24     A = ReadHBSparse("arc130.rua");
25     b = rand(130, 500);
26     t1 = 0;
27     for i = 1:nb_rep
28         tic(); x = umfpack(A,"\",b); t1 = t1 + toc();
29     end
30     t1 = t1 / nb_rep;
31     norm1 = norm(A*x-b);
32
33     t2 = 0;
34     for i = 1:nb_rep
35         tic(); x = A\b; t2 = t2 + toc();
36     end
37     t2 = t2 / nb_rep;
38     norm2 = norm(A*x-b);
39
40     disp("------------------------");
41     disp(_("Test 1: Jacobian of a nonlinear system of ordinary differential equations (ODEs) modeling a laser"));
42     disp(_("  Mean time and accuracy for umfpack (t1 and ||A*x-b||):"));
43     disp([t1 norm1]);
44     disp(_("  Time and accuracy for sparse module (t2 and ||A*x-b||):"));
45     disp([t2 norm2]);
46     disp(_("  Time comparison (t2/t1):"));
47     disp(t2/t1);
48     halt(_("\nPress Return to continue...\n"));
49
50
51     ////////////////////////////////////////////////////////////////////////////////////////////////////
52     // Test 2: Generalized eigenvalue problem Kx = (lambda)Mx
53     clear
54     nb_rep = 5;
55     B = ReadHBSparse("bcsstk24.rsa");
56     b = rand(3562, 1);
57     t1 = 0;
58     for i = 1:nb_rep
59         tic(); x = umfpack(B,"\",b); t1 = t1 + toc();
60     end
61     t1 = t1 / nb_rep;
62     norm1 = norm(B*x-b);
63
64     t2 = 0;
65     for i = 1:nb_rep
66         tic(); x = B\b; t2 = t2 + toc();
67     end
68     t2 = t2 / nb_rep;
69     norm2 = norm(B*x-b);
70
71     disp("------------------------");
72     disp(_("Test 2: Generalized eigenvalue problem Kx = (lambda)Mx"));
73     disp(_("  Mean time and accuracy for umfpack (t1 and ||A*x-b||):"));
74     disp([t1 norm1]);
75     disp(_("  Time and accuracy for sparse module (t2 and ||A*x-b||):"));
76     disp([t2 norm2]);
77     disp(_("  Time comparison (t2/t1):"));
78     disp(t2/t1);
79     halt(_("\nPress Return to continue...\n"));
80
81
82     ////////////////////////////////////////////////////////////////////////////////////////////////////
83     // Test 3: Matrices generated by the FIDAP Package
84     clear
85     nb_rep = 2;
86     C = ReadHBSparse("ex14.rua");
87     b = rand(3251, 1);
88     t1 = 0;
89     for i = 1:nb_rep
90         tic(); x = umfpack(C,"\",b); t1 = t1 + toc();
91     end
92     t1 = t1 / nb_rep;
93     norm1 = norm(C*x-b);
94
95     t2 = 0;
96     for i = 1:nb_rep
97         tic(); x = C\b; t2 = t2 + toc();
98     end
99     t2 = t2 / nb_rep;
100     norm2 = norm(C*x-b);
101
102     disp("------------------------");
103     disp(_("Test 3: Matrices generated by the FIDAP Package"));
104     disp(_("  Mean time and accuracy for umfpack (t1 and ||A*x-b||):"));
105     disp([t1 norm1]);
106     disp(_("  Time and accuracy for sparse module (t2 and ||A*x-b||):"));
107     disp([t2 norm2]);
108     disp(_("  Time comparison (t2/t1):"));
109     disp(t2/t1);
110     halt(_("\nPress Return to continue...\n"));
111
112
113     ////////////////////////////////////////////////////////////////////////////////////////////////////
114     // Test 4: Tokamak Matrices
115     clear
116     nb_rep = 1;
117     warning("off");
118     D = ReadHBSparse("utm300.rua");
119     b = rand(300, 1);
120     t1 = 0;
121     for i = 1:nb_rep
122         tic(); x = umfpack(D,"\",b); t1 = t1 + toc();
123     end
124     t1 = t1 / nb_rep;
125     norm1 = norm(D*x-b);
126
127     t2 = 0;
128     for i = 1:nb_rep
129         tic(); x = D\b; t2 = t2 + toc();
130     end
131     t2 = t2 / nb_rep;
132     norm2 = norm(D*x-b);
133
134     disp("------------------------");
135     disp(_("Test 4: Tokamak Matrices"));
136     disp(_("  Mean time and accuracy for umfpack (t1 and ||A*x-b||):"));
137     disp([t1 norm1]);
138     disp(_("  Time and accuracy for sparse module (t2 and ||A*x-b||):"));
139     disp([t2 norm2]);
140     disp(_("  Time comparison (t2/t1):"));
141     disp(t2/t1);
142     halt(_("\nPress Return to continue...\n"));
143
144
145     ////////////////////////////////////////////////////////////////////////////////////////////////////
146     // Test 5: Acoustic scattering
147     clear
148     nb_rep = 10;
149     E = ReadHBSparse("young1c.csa");
150     b = rand(841, 1);
151     t1 = 0;
152     for i = 1:nb_rep
153         tic(); x = umfpack(E,"\",b); t1 = t1 + toc();
154     end
155     t1 = t1 / nb_rep;
156     norm1 = norm(E*x-b);
157
158     t2 = 0;
159     for i = 1:nb_rep
160         tic(); x = E\b; t2 = t2 + toc();
161     end
162     t2 = t2 / nb_rep;
163     norm2 = norm(E*x-b);
164
165     disp("------------------------");
166     disp(_("Test 5: Acoustic scattering"));
167     disp(_("  Mean time and accuracy for umfpack (t1 and ||A*x-b||):"));
168     disp([t1 norm1]);
169     disp(_("  Time and accuracy for sparse module (t2 and ||A*x-b||):"));
170     disp([t2 norm2]);
171     disp(_("  Time comparison (t2/t1):"));
172     disp(t2/t1);
173     disp("------------------------");
174
175 endfunction
176
177 demos_harwell_boeing()
178 clear demos_harwell_boeing;