Annotation of OpenXM_contrib/PHC/Ada/Math_Lib/Matrices/ts_cmpmat.adb, Revision 1.1.1.1
1.1 maekawa 1: with text_io,integer_io; use text_io,integer_io;
2: with Standard_Natural_Vectors;
3: with Standard_Floating_Numbers; use Standard_Floating_Numbers;
4: with Standard_Floating_Numbers_io; use Standard_Floating_Numbers_io;
5: with Standard_Complex_Numbers;
6: with Standard_Complex_Vectors;
7: with Standard_Complex_Vectors_io;
8: with Standard_Complex_Matrices;
9: with Standard_Complex_Matrices_io;
10: with Standard_Complex_VecMats;
11: with Standard_Complex_VecMats_io;
12: with Standard_Complex_Linear_Solvers; use Standard_Complex_Linear_Solvers;
13: with Standard_Complex_Norms_Equals; use Standard_Complex_Norms_Equals;
14: with Standard_Random_Vectors; use Standard_Random_Vectors;
15: with Multprec_Random_Vectors; use Multprec_Random_Vectors;
16: with Multprec_Floating_Numbers; use Multprec_Floating_Numbers;
17: with Multprec_Floating_Numbers_io; use Multprec_Floating_Numbers_io;
18: with Multprec_Complex_Numbers;
19: with Multprec_Complex_Vectors;
20: with Multprec_Complex_Vectors_io;
21: with Multprec_Complex_Matrices;
22: with Multprec_Complex_Matrices_io;
23: with Multprec_Complex_Linear_Solvers; use Multprec_Complex_Linear_Solvers;
24: with Multprec_Complex_Norms_Equals; use Multprec_Complex_Norms_Equals;
25:
26: procedure ts_cmpmat is
27:
28: -- DESCRIPTION :
29: -- Tests the matrix packages of standard and multi-precision complex numbers.
30:
31: procedure Test_Standard_io is
32:
33: use Standard_Complex_Matrices,Standard_Complex_Matrices_io;
34:
35: n,m : natural;
36:
37: begin
38: put("Give the number of rows : "); get(n);
39: put("Give the number of columns : "); get(m);
40: declare
41: mat : Matrix(1..n,1..m);
42: begin
43: put("Give "); put(n,1); put("x"); put(m,1);
44: put_line(" complex matrix : "); get(mat);
45: put_line("Your matrix : "); put(mat); new_line;
46: end;
47: end Test_Standard_io;
48:
49: procedure Test_Standard_VecMat_io is
50:
51: use Standard_Complex_Matrices,Standard_Complex_Matrices_io;
52: use Standard_Complex_VecMats,Standard_Complex_VecMats_io;
53:
54: n,n1,n2 : natural;
55: lv : Link_to_VecMat;
56:
57: begin
58: put("Give the number of matrices : "); get(n);
59: put("Give #rows : "); get(n1);
60: put("Give #columns : "); get(n2);
61: put("Give "); put(n,1); put(" "); put(n1,1); put("-by-"); put(n2,1);
62: put_line(" integer matrices : ");
63: get(n,n1,n2,lv);
64: put_line("The vector of matrices :"); put(lv);
65: end Test_Standard_VecMat_io;
66:
67: procedure Test_Multprec_io is
68:
69: use Multprec_Complex_Matrices,Multprec_Complex_Matrices_io;
70:
71: n,m : natural;
72:
73: begin
74: put("Give the number of rows : "); get(n);
75: put("Give the number of columns : "); get(m);
76: declare
77: mat : Matrix(1..n,1..m);
78: begin
79: put("Give "); put(n,1); put("x"); put(m,1);
80: put_line(" complex matrix : "); get(mat);
81: put_line("Your matrix : "); put(mat); new_line;
82: end;
83: end Test_Multprec_io;
84:
85: function Vdm_Matrix ( v : Standard_Complex_Vectors.Vector )
86: return Standard_Complex_Matrices.Matrix is
87:
88: use Standard_Complex_Numbers;
89: use Standard_Complex_Matrices;
90:
91: n : constant natural := v'length;
92: res : Matrix(1..n,1..n);
93:
94: begin
95: for i in res'range(1) loop
96: for j in res'range(2) loop
97: res(i,j) := v(i)**(j-1);
98: end loop;
99: end loop;
100: return res;
101: end Vdm_Matrix;
102:
103: procedure lufac_Solve ( n : in natural;
104: mat : in Standard_Complex_Matrices.Matrix;
105: rhs : in Standard_Complex_Vectors.Vector ) is
106:
107: use Standard_Complex_Vectors;
108: use Standard_Complex_Vectors_io;
109: use Standard_Complex_Matrices;
110: use Standard_Complex_Matrices_io;
111:
112: wrk : Matrix(mat'range(1),mat'range(2));
113: piv : Standard_Natural_Vectors.Vector(mat'range(2));
114: info : natural;
115: res,sol,acc : Vector(rhs'range);
116: nrm : double_float;
117:
118: begin
119: put_line("Solving the linear system with lufac.");
120: wrk := mat;
121: lufac(wrk,n,piv,info);
122: put("info : "); put(info,1); new_line;
123: sol := rhs;
124: lusolve(wrk,n,piv,sol);
125: put_line("The solution vector :"); put(sol); new_line;
126: acc := mat*sol;
127: res := rhs - acc;
128: put_line("The residual : "); put(res); new_line;
129: nrm := Max_Norm(res);
130: put("Max norm of residual : "); put(nrm); new_line;
131: nrm := Sum_Norm(res);
132: put("Sum norm of residual : "); put(nrm); new_line;
133: end lufac_Solve;
134:
135: procedure lufco_Solve ( n : in natural;
136: mat : in Standard_Complex_Matrices.Matrix;
137: rhs : in Standard_Complex_Vectors.Vector ) is
138:
139: use Standard_Complex_Vectors;
140: use Standard_Complex_Vectors_io;
141: use Standard_Complex_Matrices;
142: use Standard_Complex_Matrices_io;
143:
144: wrk : Matrix(mat'range(1),mat'range(2)) := mat;
145: piv : Standard_Natural_Vectors.Vector(mat'range(2));
146: rcond,nrm : double_float;
147: res,sol : Vector(rhs'range);
148:
149: begin
150: put_line("Solving the linear system with lufco.");
151: lufco(wrk,n,piv,rcond);
152: put("inverse condition : "); put(rcond); new_line;
153: sol := rhs;
154: lusolve(wrk,n,piv,sol);
155: put_line("The solution vector :"); put(sol); new_line;
156: res := rhs - mat*sol;
157: put_line("The residual : "); put(res); new_line;
158: nrm := Max_Norm(res);
159: put("Max norm of residual : "); put(nrm); new_line;
160: nrm := Sum_Norm(res);
161: put("Sum norm of residual : "); put(nrm); new_line;
162: end lufco_Solve;
163:
164: procedure Random_Test_Standard_Linear_Solvers is
165:
166: use Standard_Complex_Vectors;
167: use Standard_Complex_Vectors_io;
168: use Standard_Complex_Matrices;
169: use Standard_Complex_Matrices_io;
170:
171: n,nb : natural;
172:
173: begin
174: new_line;
175: put_line("Testing of solving random standard linear systems.");
176: new_line;
177: put("Give the dimension : "); get(n);
178: put("Give the number of tests : "); get(nb);
179: for i in 1..nb loop
180: declare
181: mat : Matrix(1..n,1..n);
182: rhs : Vector(1..n);
183: begin
184: mat := Vdm_Matrix(Random_Vector(1,n));
185: rhs := Random_Vector(1,n);
186: lufac_Solve(n,mat,rhs);
187: lufco_Solve(n,mat,rhs);
188: end;
189: end loop;
190: end Random_Test_Standard_Linear_Solvers;
191:
192: function Vdm_Matrix ( v : Multprec_Complex_Vectors.Vector )
193: return Multprec_Complex_Matrices.Matrix is
194:
195: use Multprec_Complex_Numbers;
196: use Multprec_Complex_Matrices;
197:
198: n : constant natural := v'length;
199: res : Matrix(1..n,1..n);
200:
201: begin
202: for i in res'range(1) loop
203: for j in res'range(2) loop
204: res(i,j) := v(i)**(j-1);
205: end loop;
206: end loop;
207: return res;
208: end Vdm_Matrix;
209:
210: procedure lufac_Solve ( n : in natural;
211: mat : in Multprec_Complex_Matrices.Matrix;
212: rhs : in Multprec_Complex_Vectors.Vector ) is
213:
214: use Multprec_Complex_Vectors;
215: use Multprec_Complex_Vectors_io;
216: use Multprec_Complex_Matrices;
217: use Multprec_Complex_Matrices_io;
218:
219: wrk : Matrix(mat'range(1),mat'range(2)) := +mat;
220: piv : Standard_Natural_Vectors.Vector(mat'range(2));
221: info : natural;
222: res,sol,acc : Vector(rhs'range);
223: nrm : Floating_Number;
224:
225: begin
226: put_line("Solving the linear system with lufac.");
227: lufac(wrk,n,piv,info);
228: put("info : "); put(info,1); new_line;
229: sol := +rhs;
230: lusolve(wrk,n,piv,sol);
231: put_line("The solution vector :"); put(sol); new_line;
232: acc := mat*sol;
233: res := rhs - acc;
234: put_line("The residual : "); put(res); new_line;
235: nrm := Max_Norm(res);
236: put("Max norm of residual : "); put(nrm); new_line;
237: Clear(nrm);
238: nrm := Sum_Norm(res);
239: put("Sum norm of residual : "); put(nrm); new_line;
240: Clear(nrm);
241: Clear(wrk); Clear(acc); Clear(res); Clear(sol);
242: end lufac_Solve;
243:
244: procedure lufco_Solve ( n : in natural;
245: mat : in Multprec_Complex_Matrices.Matrix;
246: rhs : in Multprec_Complex_Vectors.Vector ) is
247:
248: use Multprec_Complex_Vectors;
249: use Multprec_Complex_Vectors_io;
250: use Multprec_Complex_Matrices;
251: use Multprec_Complex_Matrices_io;
252:
253: wrk : Matrix(mat'range(1),mat'range(2)) := +mat;
254: piv : Standard_Natural_Vectors.Vector(mat'range(2));
255: rcond,nrm : Floating_Number;
256: res,sol,acc : Vector(rhs'range);
257:
258: begin
259: put_line("Solving the linear system with lufco.");
260: lufco(wrk,n,piv,rcond);
261: put("inverse condition : "); put(rcond); new_line;
262: sol := +rhs;
263: lusolve(wrk,n,piv,sol);
264: put_line("The solution vector :"); put(sol); new_line;
265: acc := mat*sol;
266: res := rhs - acc;
267: put_line("The residual : "); put(res); new_line;
268: nrm := Max_Norm(res);
269: put("Max norm of residual : "); put(nrm); new_line;
270: Clear(nrm);
271: nrm := Sum_Norm(res);
272: put("Sum norm of residual : "); put(nrm); new_line;
273: Clear(wrk); Clear(acc); Clear(res); Clear(sol);
274: Clear(nrm); Clear(rcond);
275: end lufco_Solve;
276:
277: procedure Random_Test_Multprec_Linear_Solvers is
278:
279: use Multprec_Complex_Vectors;
280: use Multprec_Complex_Vectors_io;
281: use Multprec_Complex_Matrices;
282: use Multprec_Complex_Matrices_io;
283:
284: n,sz,nb : natural;
285:
286: begin
287: new_line;
288: put_line("Testing of solving random multi-precision linear systems.");
289: new_line;
290: put("Give the dimension : "); get(n);
291: put("Give the size of the numbers : "); get(sz);
292: put("Give the number of tests : "); get(nb);
293: for i in 1..nb loop
294: declare
295: rnd : Vector(1..n) := Random_Vector(1,n,sz);
296: mat : Matrix(1..n,1..n) := Vdm_Matrix(rnd);
297: rhs : Vector(1..n) := Random_Vector(1,n,sz);
298: begin
299: lufac_Solve(n,mat,rhs);
300: -- lufco_Solve(n,mat,rhs);
301: Clear(rnd); Clear(mat); Clear(rhs);
302: end;
303: end loop;
304: end Random_Test_Multprec_Linear_Solvers;
305:
306: procedure Main is
307:
308: ans : character;
309:
310: begin
311: new_line;
312: put_line("Interactive testing of matrices of complex numbers");
313: loop
314: new_line;
315: put_line("Choose one of the following : ");
316: put_line(" 0. exit this program. ");
317: put_line(" 1. io of matrices of standard complex numbers. ");
318: put_line(" 2. io of vectors of matrices of standard complex numbers. ");
319: put_line(" 3. test on solving random standard linear systems. ");
320: put_line(" 4. io of matrices of multi-precision complex numbers. ");
321: put_line(" 5. test on solving random multi-precision complex numbers.");
322: put("Make your choice (0,1,2,3,4 or 5) : "); get(ans);
323: exit when (ans = '0');
324: case ans is
325: when '1' => Test_Standard_io;
326: when '2' => Test_Standard_VecMat_io;
327: when '3' => Random_Test_Standard_Linear_Solvers;
328: when '4' => Test_Multprec_io;
329: when '5' => Random_Test_Multprec_Linear_Solvers;
330: when others => null;
331: end case;
332: end loop;
333: end Main;
334:
335: begin
336: Main;
337: end ts_cmpmat;
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>