![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to generic_matrices.adb CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib / PHC / Ada / Math_Lib / Matrices|
Return to generic_matrices.adb CVS log | Up to [local] / OpenXM_contrib / PHC / Ada / Math_Lib / Matrices |
1.1 maekawa 1: with unchecked_deallocation;
2:
3: package body Generic_Matrices is
4:
5: -- COMPARISON AND COPYING :
6:
7: function Equal ( a,b : Matrix ) return boolean is
8: begin
9: for i in a'range(1) loop
10: for j in a'range(2) loop
11: if not Equal(a(i,j),b(i,j))
12: then return false;
13: end if;
14: end loop;
15: end loop;
16: return true;
17: exception
18: when CONSTRAINT_ERROR => return false;
19: end Equal;
20:
21: procedure Copy ( a : in Matrix; b : in out Matrix ) is
22: begin
23: for i in a'range(1) loop
24: for j in a'range(2) loop
25: Copy(a(i,j),b(i,j));
26: end loop;
27: end loop;
28: end Copy;
29:
30: -- MATRIX-MATRIX OPERATIONS :
31:
32: function "+" ( a,b : Matrix ) return Matrix is
33:
34: res : Matrix(a'range(1),a'range(2));
35:
36: begin
37: for i in res'range(1) loop
38: for j in res'range(2) loop
39: res(i,j) := a(i,j) + b(i,j);
40: end loop;
41: end loop;
42: return res;
43: end "+";
44:
45: function "-" ( a,b : Matrix ) return Matrix is
46:
47: res : Matrix(a'range(1),a'range(2));
48:
49: begin
50: for i in res'range(1) loop
51: for j in res'range(2) loop
52: res(i,j) := a(i,j) - b(i,j);
53: end loop;
54: end loop;
55: return res;
56: end "-";
57:
58: function "+" ( a : Matrix ) return Matrix is
59:
60: res : Matrix(a'range(1),a'range(2));
61:
62: begin
63: for i in res'range(1) loop
64: for j in res'range(2) loop
65: res(i,j) := +a(i,j);
66: end loop;
67: end loop;
68: return res;
69: end "+";
70:
71: function "-" ( a : Matrix ) return Matrix is
72:
73: res : Matrix(a'range(1),a'range(2));
74:
75: begin
76: for i in res'range(1) loop
77: for j in res'range(2) loop
78: res(i,j) := -a(i,j);
79: end loop;
80: end loop;
81: return res;
82: end "-";
83:
84: function "*" ( a,b : Matrix ) return Matrix is
85:
86: res : Matrix(a'range(1),b'range(2));
87: acc : number;
88:
89: begin
90: for i in res'range(1) loop
91: for j in res'range(2) loop
92: Copy(a(i,a'first(2))*b(b'first(1),j),res(i,j));
93: for k in a'first(2)+1..a'last(2) loop
94: acc := a(i,k)*b(k,j);
95: Add(res(i,j),acc);
96: Clear(acc);
97: end loop;
98: end loop;
99: end loop;
100: return res;
101: end "*";
102:
103: procedure Mul1 ( a : in out Matrix; b : in Matrix ) is
104:
105: temp : Vector(a'range(2));
106: acc : number;
107:
108: begin
109: for i in a'range(1) loop
110: for j in b'range(2) loop
111: Copy(a(i,a'first(2))*b(b'first(1),j),temp(j));
112: for k in a'first(2)+1..a'last(2) loop
113: acc := a(i,k)*b(k,j);
114: Add(temp(j),acc);
115: Clear(acc);
116: end loop;
117: end loop;
118: for j in a'range(2) loop
119: Copy(temp(j),a(i,j));
120: end loop;
121: end loop;
122: end Mul1;
123:
124: procedure Mul2 ( a : in Matrix; b : in out Matrix ) is
125:
126: temp : Vector(a'range(1));
127: acc : number;
128:
129: begin
130: for i in b'range(2) loop
131: for j in a'range(1) loop
132: Copy(a(j,a'first(1))*b(a'first(1),i),temp(j));
133: for k in a'first(1)+1..a'last(1) loop
134: acc := a(j,k)*b(k,i);
135: Add(temp(j),acc);
136: Clear(acc);
137: end loop;
138: end loop;
139: for j in b'range(1) loop
140: Copy(temp(j),b(j,i));
141: end loop;
142: end loop;
143: end Mul2;
144:
145: -- MATRIX-VECTOR OPERATIONS :
146:
147: function "*" ( a : Matrix; v : Vector ) return Vector is
148:
149: res : Vector(a'range(1));
150: acc : number;
151:
152: begin
153: for i in res'range loop
154: Copy(a(i,a'first(2))*v(v'first),res(i));
155: for j in a'first(2)+1..a'last(2) loop
156: acc := a(i,j)*v(j);
157: Add(res(i),acc);
158: Clear(acc);
159: end loop;
160: end loop;
161: return res;
162: end "*";
163:
164: function "*" ( v : Vector; a : Matrix ) return Vector is
165:
166: res : Vector(a'range(2));
167: acc : number;
168:
169: begin
170: for j in res'range loop
171: Copy(v(v'first)*a(a'first(1),j),res(j));
172: for i in a'first(1)+1..a'last(1) loop
173: acc := v(i)*a(i,j);
174: Add(res(j),acc);
175: Clear(acc);
176: end loop;
177: end loop;
178: return res;
179: end "*";
180:
181: procedure Mul ( a : in Matrix; v : in out Vector ) is
182:
183: tv : Vector(v'range);
184: acc : number;
185:
186: begin
187: for i in v'range loop
188: Copy(a(i,a'first(2))*v(v'first),tv(i));
189: for j in a'first(2)+1..a'last(2) loop
190: acc := a(i,j)*v(j);
191: Add(tv(i),a(i,j)*v(j));
192: Clear(acc);
193: end loop;
194: end loop;
195: for i in v'range loop
196: v(i) := tv(i);
197: end loop;
198: end Mul;
199:
200: procedure Mul ( v : in out Vector; a : in Matrix ) is
201:
202: tv : Vector(v'range);
203: acc : number;
204:
205: begin
206: for j in v'range loop
207: Copy(v(v'first)*a(a'first(1),j),tv(j));
208: for i in a'first(1)+1..a'last(1) loop
209: acc := v(j)*a(i,j);
210: Add(tv(j),acc);
211: Clear(acc);
212: end loop;
213: end loop;
214: for i in v'range loop
215: v(i) := tv(i);
216: end loop;
217: end Mul;
218:
219: -- DESTRUCTORS :
220:
221: procedure Clear ( a : in out Matrix ) is
222: begin
223: for i in a'range(1) loop
224: for j in a'range(2) loop
225: Clear(a(i,j));
226: end loop;
227: end loop;
228: end Clear;
229:
230: procedure Clear ( a : in out Link_to_Matrix ) is
231:
232: procedure free is new unchecked_deallocation(Matrix,Link_to_Matrix);
233:
234: begin
235: if a /= null
236: then Clear(a.all);
237: end if;
238: free(a);
239: end Clear;
240:
241: end Generic_Matrices;