1 ------------------------------------------------------------------------------
3 -- GNAT RUN-TIME COMPONENTS --
5 -- S Y S T E M . G E N E R I C _ A R R A Y _ O P E R A T I O N S --
9 -- Copyright (C) 2006-2009, Free Software Foundation, Inc. --
11 -- GNAT is free software; you can redistribute it and/or modify it under --
12 -- terms of the GNU General Public License as published by the Free Soft- --
13 -- ware Foundation; either version 3, or (at your option) any later ver- --
14 -- sion. GNAT is distributed in the hope that it will be useful, but WITH- --
15 -- OUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY --
16 -- or FITNESS FOR A PARTICULAR PURPOSE. --
18 -- As a special exception under Section 7 of GPL version 3, you are granted --
19 -- additional permissions described in the GCC Runtime Library Exception, --
20 -- version 3.1, as published by the Free Software Foundation. --
22 -- You should have received a copy of the GNU General Public License and --
23 -- a copy of the GCC Runtime Library Exception along with this program; --
24 -- see the files COPYING3 and COPYING.RUNTIME respectively. If not, see --
25 -- <http://www.gnu.org/licenses/>. --
27 -- GNAT was originally developed by the GNAT team at New York University. --
28 -- Extensive contributions were provided by Ada Core Technologies Inc. --
30 ------------------------------------------------------------------------------
32 package body System.Generic_Array_Operations is
34 -- The local function Check_Unit_Last computes the index
35 -- of the last element returned by Unit_Vector or Unit_Matrix.
36 -- A separate function is needed to allow raising Constraint_Error
37 -- before declaring the function result variable. The result variable
38 -- needs to be declared first, to allow front-end inlining.
40 function Check_Unit_Last
43 First : Integer) return Integer;
44 pragma Inline_Always (Check_Unit_Last);
46 function Square_Matrix_Length (A : Matrix) return Natural is
48 if A'Length (1) /= A'Length (2) then
49 raise Constraint_Error with "matrix is not square";
53 end Square_Matrix_Length;
59 function Check_Unit_Last
62 First : Integer) return Integer is
64 -- Order the tests carefully to avoid overflow
67 or else First > Integer'Last - Order + 1
68 or else Index > First + (Order - 1)
70 raise Constraint_Error;
73 return First + (Order - 1);
80 function Inner_Product
85 R : Result_Scalar := Zero;
88 if Left'Length /= Right'Length then
89 raise Constraint_Error with
90 "vectors are of different length in inner product";
93 for J in Left'Range loop
94 R := R + Left (J) * Right (J - Left'First + Right'First);
100 ----------------------------------
101 -- Matrix_Elementwise_Operation --
102 ----------------------------------
104 function Matrix_Elementwise_Operation (X : X_Matrix) return Result_Matrix is
105 R : Result_Matrix (X'Range (1), X'Range (2));
108 for J in R'Range (1) loop
109 for K in R'Range (2) loop
110 R (J, K) := Operation (X (J, K));
115 end Matrix_Elementwise_Operation;
117 ----------------------------------
118 -- Vector_Elementwise_Operation --
119 ----------------------------------
121 function Vector_Elementwise_Operation (X : X_Vector) return Result_Vector is
122 R : Result_Vector (X'Range);
125 for J in R'Range loop
126 R (J) := Operation (X (J));
130 end Vector_Elementwise_Operation;
132 -----------------------------------------
133 -- Matrix_Matrix_Elementwise_Operation --
134 -----------------------------------------
136 function Matrix_Matrix_Elementwise_Operation
138 Right : Right_Matrix)
141 R : Result_Matrix (Left'Range (1), Left'Range (2));
143 if Left'Length (1) /= Right'Length (1)
144 or else Left'Length (2) /= Right'Length (2)
146 raise Constraint_Error with
147 "matrices are of different dimension in elementwise operation";
150 for J in R'Range (1) loop
151 for K in R'Range (2) loop
156 (J - R'First (1) + Right'First (1),
157 K - R'First (2) + Right'First (2)));
162 end Matrix_Matrix_Elementwise_Operation;
164 ------------------------------------------------
165 -- Matrix_Matrix_Scalar_Elementwise_Operation --
166 ------------------------------------------------
168 function Matrix_Matrix_Scalar_Elementwise_Operation
171 Z : Z_Scalar) return Result_Matrix
173 R : Result_Matrix (X'Range (1), X'Range (2));
176 if X'Length (1) /= Y'Length (1)
177 or else X'Length (2) /= Y'Length (2)
179 raise Constraint_Error with
180 "matrices are of different dimension in elementwise operation";
183 for J in R'Range (1) loop
184 for K in R'Range (2) loop
188 Y (J - R'First (1) + Y'First (1),
189 K - R'First (2) + Y'First (2)),
195 end Matrix_Matrix_Scalar_Elementwise_Operation;
197 -----------------------------------------
198 -- Vector_Vector_Elementwise_Operation --
199 -----------------------------------------
201 function Vector_Vector_Elementwise_Operation
203 Right : Right_Vector) return Result_Vector
205 R : Result_Vector (Left'Range);
208 if Left'Length /= Right'Length then
209 raise Constraint_Error with
210 "vectors are of different length in elementwise operation";
213 for J in R'Range loop
214 R (J) := Operation (Left (J), Right (J - R'First + Right'First));
218 end Vector_Vector_Elementwise_Operation;
220 ------------------------------------------------
221 -- Vector_Vector_Scalar_Elementwise_Operation --
222 ------------------------------------------------
224 function Vector_Vector_Scalar_Elementwise_Operation
227 Z : Z_Scalar) return Result_Vector
229 R : Result_Vector (X'Range);
232 if X'Length /= Y'Length then
233 raise Constraint_Error with
234 "vectors are of different length in elementwise operation";
237 for J in R'Range loop
238 R (J) := Operation (X (J), Y (J - X'First + Y'First), Z);
242 end Vector_Vector_Scalar_Elementwise_Operation;
244 -----------------------------------------
245 -- Matrix_Scalar_Elementwise_Operation --
246 -----------------------------------------
248 function Matrix_Scalar_Elementwise_Operation
250 Right : Right_Scalar) return Result_Matrix
252 R : Result_Matrix (Left'Range (1), Left'Range (2));
255 for J in R'Range (1) loop
256 for K in R'Range (2) loop
257 R (J, K) := Operation (Left (J, K), Right);
262 end Matrix_Scalar_Elementwise_Operation;
264 -----------------------------------------
265 -- Vector_Scalar_Elementwise_Operation --
266 -----------------------------------------
268 function Vector_Scalar_Elementwise_Operation
270 Right : Right_Scalar) return Result_Vector
272 R : Result_Vector (Left'Range);
275 for J in R'Range loop
276 R (J) := Operation (Left (J), Right);
280 end Vector_Scalar_Elementwise_Operation;
282 -----------------------------------------
283 -- Scalar_Matrix_Elementwise_Operation --
284 -----------------------------------------
286 function Scalar_Matrix_Elementwise_Operation
288 Right : Right_Matrix) return Result_Matrix
290 R : Result_Matrix (Right'Range (1), Right'Range (2));
293 for J in R'Range (1) loop
294 for K in R'Range (2) loop
295 R (J, K) := Operation (Left, Right (J, K));
300 end Scalar_Matrix_Elementwise_Operation;
302 -----------------------------------------
303 -- Scalar_Vector_Elementwise_Operation --
304 -----------------------------------------
306 function Scalar_Vector_Elementwise_Operation
308 Right : Right_Vector) return Result_Vector
310 R : Result_Vector (Right'Range);
313 for J in R'Range loop
314 R (J) := Operation (Left, Right (J));
318 end Scalar_Vector_Elementwise_Operation;
320 ---------------------------
321 -- Matrix_Matrix_Product --
322 ---------------------------
324 function Matrix_Matrix_Product
326 Right : Right_Matrix) return Result_Matrix
328 R : Result_Matrix (Left'Range (1), Right'Range (2));
331 if Left'Length (2) /= Right'Length (1) then
332 raise Constraint_Error with
333 "incompatible dimensions in matrix multiplication";
336 for J in R'Range (1) loop
337 for K in R'Range (2) loop
339 S : Result_Scalar := Zero;
341 for M in Left'Range (2) loop
343 * Right (M - Left'First (2) + Right'First (1), K);
352 end Matrix_Matrix_Product;
354 ---------------------------
355 -- Matrix_Vector_Product --
356 ---------------------------
358 function Matrix_Vector_Product
360 Right : Right_Vector) return Result_Vector
362 R : Result_Vector (Left'Range (1));
365 if Left'Length (2) /= Right'Length then
366 raise Constraint_Error with
367 "incompatible dimensions in matrix-vector multiplication";
370 for J in Left'Range (1) loop
372 S : Result_Scalar := Zero;
374 for K in Left'Range (2) loop
375 S := S + Left (J, K) * Right (K - Left'First (2) + Right'First);
383 end Matrix_Vector_Product;
389 function Outer_Product
391 Right : Right_Vector) return Matrix
393 R : Matrix (Left'Range, Right'Range);
396 for J in R'Range (1) loop
397 for K in R'Range (2) loop
398 R (J, K) := Left (J) * Right (K);
409 procedure Transpose (A : Matrix; R : out Matrix) is
411 for J in R'Range (1) loop
412 for K in R'Range (2) loop
413 R (J, K) := A (K - R'First (2) + A'First (1),
414 J - R'First (1) + A'First (2));
419 -------------------------------
420 -- Update_Matrix_With_Matrix --
421 -------------------------------
423 procedure Update_Matrix_With_Matrix (X : in out X_Matrix; Y : Y_Matrix) is
425 if X'Length (1) /= Y'Length (1)
426 or else X'Length (2) /= Y'Length (2)
428 raise Constraint_Error with
429 "matrices are of different dimension in update operation";
432 for J in X'Range (1) loop
433 for K in X'Range (2) loop
434 Update (X (J, K), Y (J - X'First (1) + Y'First (1),
435 K - X'First (2) + Y'First (2)));
438 end Update_Matrix_With_Matrix;
440 -------------------------------
441 -- Update_Vector_With_Vector --
442 -------------------------------
444 procedure Update_Vector_With_Vector (X : in out X_Vector; Y : Y_Vector) is
446 if X'Length /= Y'Length then
447 raise Constraint_Error with
448 "vectors are of different length in update operation";
451 for J in X'Range loop
452 Update (X (J), Y (J - X'First + Y'First));
454 end Update_Vector_With_Vector;
462 First_1 : Integer := 1;
463 First_2 : Integer := 1) return Matrix
465 R : Matrix (First_1 .. Check_Unit_Last (First_1, Order, First_1),
466 First_2 .. Check_Unit_Last (First_2, Order, First_2));
469 R := (others => (others => Zero));
471 for J in 0 .. Order - 1 loop
472 R (First_1 + J, First_2 + J) := One;
485 First : Integer := 1) return Vector
487 R : Vector (First .. Check_Unit_Last (Index, Order, First));
489 R := (others => Zero);
494 ---------------------------
495 -- Vector_Matrix_Product --
496 ---------------------------
498 function Vector_Matrix_Product
500 Right : Matrix) return Result_Vector
502 R : Result_Vector (Right'Range (2));
505 if Left'Length /= Right'Length (2) then
506 raise Constraint_Error with
507 "incompatible dimensions in vector-matrix multiplication";
510 for J in Right'Range (2) loop
512 S : Result_Scalar := Zero;
515 for K in Right'Range (1) loop
516 S := S + Left (J - Right'First (1) + Left'First) * Right (K, J);
524 end Vector_Matrix_Product;
526 end System.Generic_Array_Operations;