Annotation of OpenXM_contrib/PHC/Ada/Math_Lib/Matrices/ts_fltmat.adb, Revision 1.1.1.1
1.1 maekawa 1: with text_io,integer_io; use text_io,integer_io;
2: with Multprec_Integer_Numbers; use Multprec_Integer_Numbers;
3: with Standard_Natural_Vectors;
4: with Standard_Floating_Numbers; use Standard_Floating_Numbers;
5: with Standard_Floating_Numbers_io; use Standard_Floating_Numbers_io;
6: --with Standard_Random_Numbers; use Standard_Random_Numbers;
7: with Standard_Floating_Vectors;
8: with Standard_Floating_Vectors_io;
9: with Standard_Floating_Matrices;
10: with Standard_Floating_Matrices_io;
11: with Standard_Floating_VecMats;
12: with Standard_Floating_VecMats_io;
13: with Standard_Floating_Linear_Solvers; use Standard_Floating_Linear_Solvers;
14: with Standard_Floating_Norms_Equals; use Standard_Floating_Norms_Equals;
15: with Standard_Random_Vectors; use Standard_Random_Vectors;
16: with Standard_Random_Matrices; use Standard_Random_Matrices;
17: with Multprec_Random_Vectors; use Multprec_Random_Vectors;
18: with Multprec_Random_Matrices; use Multprec_Random_Matrices;
19: with Multprec_Floating_Numbers; use Multprec_Floating_Numbers;
20: with Multprec_Floating_Numbers_io; use Multprec_Floating_Numbers_io;
21: --with Multprec_Random_Numbers; use Multprec_Random_Numbers;
22: with Multprec_Floating_Vectors;
23: with Multprec_Floating_Vectors_io;
24: with Multprec_Floating_Matrices;
25: with Multprec_Floating_Matrices_io;
26: with Multprec_Floating_Linear_Solvers; use Multprec_Floating_Linear_Solvers;
27: with Multprec_Floating_Norms_Equals; use Multprec_Floating_Norms_Equals;
28:
29: procedure ts_fltmat is
30:
31: -- DESCRIPTION :
32: -- Tests the matrix packages of standard and multi-precision floats.
33:
34: function Vdm_Matrix ( v : Multprec_Floating_Vectors.Vector )
35: return Multprec_Floating_Matrices.Matrix is
36:
37: use Multprec_Floating_Matrices;
38:
39: n : constant natural := v'length;
40: res : Matrix(1..n,1..n);
41:
42: begin
43: for i in res'range(1) loop
44: for j in res'range(2) loop
45: res(i,j) := v(i)**(j-1);
46: end loop;
47: end loop;
48: return res;
49: end Vdm_Matrix;
50:
51: function Vdm_Matrix ( v : Standard_Floating_Vectors.Vector )
52: return Standard_Floating_Matrices.Matrix is
53:
54: use Standard_Floating_Matrices;
55:
56: n : constant natural := v'length;
57: res : Matrix(1..n,1..n);
58:
59: begin
60: for i in res'range(1) loop
61: for j in res'range(2) loop
62: res(i,j) := v(i)**(j-1);
63: end loop;
64: end loop;
65: return res;
66: end Vdm_Matrix;
67:
68: procedure Write ( f : in Multprec_Floating_Numbers.Floating_Number ) is
69:
70: frac : Integer_Number := Fraction(f);
71: expo : Integer_Number := Exponent(f);
72:
73: begin
74: put("Fraction :");
75: for i in reverse 0..Size(frac) loop
76: put(Coefficient(frac,i));
77: end loop;
78: new_line;
79: put("Exponent :");
80: for i in reverse 0..Size(expo) loop
81: put(Coefficient(expo,i));
82: end loop;
83: new_line;
84: end Write;
85:
86: procedure Write ( v : in Multprec_Floating_Vectors.Vector ) is
87: begin
88: for i in v'range loop
89: Write(v(i)); new_line;
90: end loop;
91: end Write;
92:
93: procedure Write ( m : in Multprec_Floating_Matrices.Matrix ) is
94: begin
95: for i in m'range(1) loop
96: for j in m'range(2) loop
97: Write(m(i,j)); new_line;
98: end loop;
99: end loop;
100: end Write;
101:
102: procedure Set_Size ( v : in out Multprec_Floating_Vectors.Vector;
103: k : in natural ) is
104:
105: -- DESCRIPTION :
106: -- Sets the size of the elements in v to k.
107:
108: use Multprec_Floating_Numbers,Multprec_Floating_Vectors;
109:
110: begin
111: for i in v'range loop
112: Set_Size(v(i),k);
113: end loop;
114: --Write(v);
115: end Set_Size;
116:
117: procedure Set_Size ( m : in out Multprec_Floating_Matrices.Matrix;
118: k : in natural ) is
119:
120: -- DESCRIPTION :
121: -- Sets the size of the elements in m to k.
122:
123: use Multprec_Floating_Numbers,Multprec_Floating_Matrices;
124:
125: begin
126: for i in m'range(1) loop
127: for j in m'range(2) loop
128: Set_Size(m(i,j),k);
129: end loop;
130: end loop;
131: --Write(m);
132: end Set_Size;
133:
134: procedure Test_Standard_io is
135:
136: use Standard_Floating_Matrices,Standard_Floating_Matrices_io;
137:
138: n,m : natural;
139:
140: begin
141: put("Give the number of rows : "); get(n);
142: put("Give the number of columns : "); get(m);
143: declare
144: mat : Matrix(1..n,1..m);
145: begin
146: put("Give "); put(n,1); put("x"); put(m,1);
147: put_line(" floating matrix : "); get(mat);
148: put_line("Your matrix : "); put(mat); new_line;
149: end;
150: end Test_Standard_io;
151:
152: procedure Test_Standard_VecMat_io is
153:
154: use Standard_Floating_Matrices,Standard_Floating_Matrices_io;
155: use Standard_Floating_VecMats,Standard_Floating_VecMats_io;
156:
157: n,n1,n2 : natural;
158: lv : Link_to_VecMat;
159:
160: begin
161: put("Give the number of matrices : "); get(n);
162: put("Give #rows : "); get(n1);
163: put("Give #columns : "); get(n2);
164: put("Give "); put(n,1); put(" "); put(n1,1); put("-by-"); put(n2,1);
165: put_line(" floating matrices : ");
166: get(n,n1,n2,lv);
167: put_line("The vector of matrices :"); put(lv);
168: end Test_Standard_VecMat_io;
169:
170: procedure lufco_Solve ( n : in natural;
171: mat : in Standard_Floating_Matrices.Matrix;
172: rhs : in Standard_Floating_Vectors.Vector ) is
173:
174: use Standard_Floating_Vectors;
175: use Standard_Floating_Vectors_io;
176: use Standard_Floating_Matrices;
177: use Standard_Floating_Matrices_io;
178:
179: wrk : Matrix(mat'range(1),mat'range(2)) := mat;
180: piv : Standard_Natural_Vectors.Vector(mat'range(2));
181: rcond,nrm : double_float;
182: res,sol : Vector(rhs'range);
183:
184: begin
185: put_line("Solving the linear system with lufco.");
186: lufco(wrk,n,piv,rcond);
187: put("inverse condition : "); put(rcond); new_line;
188: sol := rhs;
189: lusolve(wrk,n,piv,sol);
190: put_line("The solution vector :"); put(sol); new_line;
191: res := rhs - mat*sol;
192: put_line("The residual : "); put(res); new_line;
193: nrm := Max_Norm(res);
194: put("Max norm of residual : "); put(nrm); new_line;
195: nrm := Sum_Norm(res);
196: put("Sum norm of residual : "); put(nrm); new_line;
197: end lufco_Solve;
198:
199: procedure lufac_Solve ( n : in natural;
200: mat : in Standard_Floating_Matrices.Matrix;
201: rhs : in Standard_Floating_Vectors.Vector ) is
202:
203: use Standard_Floating_Vectors;
204: use Standard_Floating_Vectors_io;
205: use Standard_Floating_Matrices;
206: use Standard_Floating_Matrices_io;
207:
208: wrk : Matrix(mat'range(1),mat'range(2)) := mat;
209: piv : Standard_Natural_Vectors.Vector(mat'range(2));
210: info : natural;
211: res,sol : Vector(rhs'range);
212: nrm : double_float;
213:
214: begin
215: put_line("Solving the linear system with lufac.");
216: lufac(wrk,n,piv,info);
217: put("info : "); put(info,1); new_line;
218: sol := rhs;
219: lusolve(wrk,n,piv,sol);
220: put_line("The solution vector :"); put(sol); new_line;
221: res := rhs - mat*sol;
222: put_line("The residual : "); put(res); new_line;
223: nrm := Max_Norm(res);
224: put("Max norm of residual : "); put(nrm); new_line;
225: nrm := Sum_Norm(res);
226: put("Sum norm of residual : "); put(nrm); new_line;
227: end lufac_Solve;
228:
229: procedure Interactive_Test_Standard_Linear_Solvers is
230:
231: use Standard_Floating_Vectors;
232: use Standard_Floating_Vectors_io;
233: use Standard_Floating_Matrices;
234: use Standard_Floating_Matrices_io;
235:
236: n : natural;
237:
238: begin
239: new_line;
240: put_line("Interactive testing of solving standard linear systems.");
241: new_line;
242: put("Give the dimension : "); get(n);
243: declare
244: mat : Matrix(1..n,1..n);
245: rhs : Vector(1..n);
246: begin
247: put("Give "); put(n,1); put("x"); put(n,1);
248: put_line(" floating matrix : "); get(mat);
249: put_line("-> the matrix : "); put(mat);
250: put("Give "); put(n,1); put_line(" floating-numbers : "); get(rhs);
251: put_line("-> right-hand side vector : "); put(rhs); new_line;
252: lufac_Solve(n,mat,rhs);
253: lufco_Solve(n,mat,rhs);
254: end;
255: end Interactive_Test_Standard_Linear_Solvers;
256:
257: procedure Random_Test_Standard_Linear_Solvers is
258:
259: use Standard_Floating_Vectors;
260: use Standard_Floating_Vectors_io;
261: use Standard_Floating_Matrices;
262: use Standard_Floating_Matrices_io;
263:
264: n,nb : natural;
265:
266: begin
267: new_line;
268: put_line("Testing of solving random standard linear systems.");
269: new_line;
270: put("Give the dimension : "); get(n);
271: put("Give the number of tests : "); get(nb);
272: for i in 1..nb loop
273: declare
274: mat : Matrix(1..n,1..n);
275: rhs : Vector(1..n);
276: begin
277: mat := Vdm_Matrix(Random_Vector(1,n));
278: rhs := Random_Vector(1,n);
279: lufac_Solve(n,mat,rhs);
280: lufco_Solve(n,mat,rhs);
281: end;
282: end loop;
283: end Random_Test_Standard_Linear_Solvers;
284:
285: procedure Test_Multprec_io is
286:
287: use Multprec_Floating_Matrices,Multprec_Floating_Matrices_io;
288:
289: n,m : natural;
290:
291: begin
292: put("Give the number of rows : "); get(n);
293: put("Give the number of columns : "); get(m);
294: declare
295: mat : Matrix(1..n,1..m);
296: begin
297: put("Give "); put(n,1); put("x"); put(m,1);
298: put_line(" floating matrix : "); get(mat);
299: put_line("Your matrix : "); put(mat); new_line;
300: end;
301: end Test_Multprec_io;
302:
303: procedure lufco_Solve ( n : in natural;
304: mat : in Multprec_Floating_Matrices.Matrix;
305: rhs : in Multprec_Floating_Vectors.Vector ) is
306:
307: use Multprec_Floating_Vectors;
308: use Multprec_Floating_Vectors_io;
309: use Multprec_Floating_Matrices;
310: use Multprec_Floating_Matrices_io;
311:
312: wrk : Matrix(mat'range(1),mat'range(2));
313: piv : Standard_Natural_Vectors.Vector(mat'range(2));
314: rcond,nrm : Floating_Number;
315: res,sol,acc : Vector(rhs'range);
316:
317: begin
318: put_line("Solving the linear system with lufco.");
319: Copy(mat,wrk);
320: lufco(wrk,n,piv,rcond);
321: put("inverse condition : "); put(rcond); new_line;
322: Copy(rhs,sol);
323: lusolve(wrk,n,piv,sol);
324: put_line("The solution vector :"); put(sol); new_line;
325: acc := mat*sol;
326: res := rhs - acc;
327: put_line("The residual : "); put(res); new_line;
328: nrm := Max_Norm(res);
329: put("Max norm of residual : "); put(nrm); new_line;
330: Clear(nrm);
331: nrm := Sum_Norm(res);
332: put("Sum norm of residual : "); put(nrm); new_line;
333: Clear(nrm); Clear(rcond);
334: Clear(wrk); Clear(acc); Clear(res); Clear(sol);
335: end lufco_Solve;
336:
337: procedure lufac_Solve ( n : in natural;
338: mat : in Multprec_Floating_Matrices.Matrix;
339: rhs : in Multprec_Floating_Vectors.Vector ) is
340:
341: use Multprec_Floating_Vectors;
342: use Multprec_Floating_Vectors_io;
343: use Multprec_Floating_Matrices;
344: use Multprec_Floating_Matrices_io;
345:
346: wrk : Matrix(mat'range(1),mat'range(2));
347: piv : Standard_Natural_Vectors.Vector(mat'range(2));
348: info : natural;
349: res,sol,acc : Vector(rhs'range);
350: nrm : Floating_Number;
351:
352: begin
353: put_line("Solving the linear system with lufac.");
354: Copy(mat,wrk);
355: lufac(wrk,n,piv,info);
356: put("info : "); put(info,1); new_line;
357: Copy(rhs,sol);
358: lusolve(wrk,n,piv,sol);
359: put_line("The solution vector :"); put(sol); new_line;
360: acc := mat*sol;
361: res := rhs - acc;
362: put_line("The residual : "); put(res); new_line;
363: nrm := Max_Norm(res);
364: put("Max norm of residual : "); put(nrm); new_line;
365: Clear(nrm);
366: nrm := Sum_Norm(res);
367: put("Sum norm of residual : "); put(nrm); new_line;
368: Clear(nrm); Clear(wrk); Clear(acc); Clear(res); Clear(sol);
369: end lufac_Solve;
370:
371: procedure Interactive_Test_Multprec_Linear_Solvers is
372:
373: use Multprec_Floating_Vectors;
374: use Multprec_Floating_Vectors_io;
375: use Multprec_Floating_Matrices;
376: use Multprec_Floating_Matrices_io;
377:
378: ans : character;
379: n : natural;
380:
381: begin
382: put("Give the dimension : "); get(n);
383: declare
384: mat : Matrix(1..n,1..n);
385: rhs : Vector(1..n);
386: sz : integer;
387: begin
388: put("Give "); put(n,1); put("x"); put(n,1);
389: put_line(" floating matrix : "); get(mat);
390: put_line("-> the matrix : "); put(mat);
391: put("Give "); put(n,1); put_line(" floating-numbers : "); get(rhs);
392: put_line("-> right-hand side vector : "); put(rhs); new_line;
393: put("Give the size (-1 for default) : "); get(sz);
394: if sz >= 0
395: then Set_Size(mat,sz);
396: Set_Size(rhs,sz);
397: end if;
398: loop
399: lufac_Solve(n,mat,rhs);
400: lufco_Solve(n,mat,rhs);
401: put("Do you want to resolve with other precision ? (y/n) "); get(ans);
402: exit when ans /= 'y';
403: put("Give the size : "); get(sz);
404: Set_Size(mat,sz);
405: Set_Size(rhs,sz);
406: end loop;
407: end;
408: end Interactive_Test_Multprec_Linear_Solvers;
409:
410: procedure Random_Test_Multprec_Linear_Solvers is
411:
412: use Multprec_Floating_Vectors;
413: use Multprec_Floating_Vectors_io;
414: use Multprec_Floating_Matrices;
415: use Multprec_Floating_Matrices_io;
416:
417: n,sz,nb : natural;
418:
419: begin
420: put("Give the dimension : "); get(n);
421: put("Give the size of the numbers : "); get(sz);
422: put("Give the number of tests : "); get(nb);
423: for i in 1..nb loop
424: declare
425: mat : Matrix(1..n,1..n);
426: rhs : Vector(1..n);
427: begin
428: mat := Vdm_Matrix(Random_Vector(1,n,sz)); --Random_Matrix(n,n,sz);
429: rhs := Random_Vector(n,sz);
430: -- lufac_Solve(n,mat,rhs);
431: lufco_Solve(n,mat,rhs);
432: Clear(mat); Clear(rhs);
433: end;
434: end loop;
435: end Random_Test_Multprec_Linear_Solvers;
436:
437: procedure Main is
438:
439: ans : character;
440:
441: begin
442: new_line;
443: put_line("Interactive testing of matrices of floating numbers");
444: loop
445: new_line;
446: put_line("Choose one of the following : ");
447: put_line(" 0. exit this program.");
448: put_line(" 1. io of matrices of standard numbers.");
449: put_line(" 2. io of vectos of matrices of standard numbers.");
450: put_line(" 3. interactive test on solving standard linear systems.");
451: put_line(" 4. test on solving random standard linear systems.");
452: put_line(" 5. io of matrices of multi-precision numbers.");
453: put_line(" 6. interactive test on solving multi-precision systems.");
454: put_line(" 7. test on solving random multi-precision systems.");
455: put("Make your choice (0,1,2,3,4,5,6 or 7) : "); get(ans);
456: exit when ans = '0';
457: case ans is
458: when '1' => Test_Standard_io;
459: when '2' => Test_Standard_VecMat_io;
460: when '3' => Interactive_Test_Standard_Linear_Solvers;
461: when '4' => Random_Test_Standard_Linear_Solvers;
462: when '5' => Test_Multprec_io;
463: when '6' => Interactive_Test_Multprec_Linear_Solvers;
464: when '7' => Random_Test_Multprec_Linear_Solvers;
465: when others => null;
466: end case;
467: end loop;
468: end Main;
469:
470: begin
471: Main;
472: end ts_fltmat;
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>