2014-02-20 Richard Biener <rguenther@suse.de>
[official-gcc.git] / gcc / ada / a-ngcoar.adb
blobca0c58c36f222808e679f3b91761c7b1fd68eb0b
1 ------------------------------------------------------------------------------
2 -- --
3 -- GNAT COMPILER COMPONENTS --
4 -- --
5 -- ADA.NUMERICS.GENERIC_COMPLEX_ARRAYS --
6 -- --
7 -- B o d y --
8 -- --
9 -- Copyright (C) 2006-2012, Free Software Foundation, Inc. --
10 -- --
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. --
17 -- --
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. --
21 -- --
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/>. --
26 -- --
27 -- GNAT was originally developed by the GNAT team at New York University. --
28 -- Extensive contributions were provided by Ada Core Technologies Inc. --
29 -- --
30 ------------------------------------------------------------------------------
32 with System.Generic_Array_Operations; use System.Generic_Array_Operations;
33 with Ada.Numerics; use Ada.Numerics;
35 package body Ada.Numerics.Generic_Complex_Arrays is
37 -- Operations that are defined in terms of operations on the type Real,
38 -- such as addition, subtraction and scaling, are computed in the canonical
39 -- way looping over all elements.
41 package Ops renames System.Generic_Array_Operations;
43 subtype Real is Real_Arrays.Real;
44 -- Work around visibility bug ???
46 function Is_Non_Zero (X : Complex) return Boolean is (X /= (0.0, 0.0));
47 -- Needed by Back_Substitute
49 procedure Back_Substitute is new Ops.Back_Substitute
50 (Scalar => Complex,
51 Matrix => Complex_Matrix,
52 Is_Non_Zero => Is_Non_Zero);
54 procedure Forward_Eliminate is new Ops.Forward_Eliminate
55 (Scalar => Complex,
56 Real => Real'Base,
57 Matrix => Complex_Matrix,
58 Zero => (0.0, 0.0),
59 One => (1.0, 0.0));
61 procedure Transpose is new Ops.Transpose
62 (Scalar => Complex,
63 Matrix => Complex_Matrix);
65 -- Helper function that raises a Constraint_Error is the argument is
66 -- not a square matrix, and otherwise returns its length.
68 function Length is new Square_Matrix_Length (Complex, Complex_Matrix);
70 -- Instant a generic square root implementation here, in order to avoid
71 -- instantiating a complete copy of Generic_Elementary_Functions.
72 -- Speed of the square root is not a big concern here.
74 function Sqrt is new Ops.Sqrt (Real'Base);
76 -- Instantiating the following subprograms directly would lead to
77 -- name clashes, so use a local package.
79 package Instantiations is
81 ---------
82 -- "*" --
83 ---------
85 function "*" is new Vector_Scalar_Elementwise_Operation
86 (Left_Scalar => Complex,
87 Right_Scalar => Complex,
88 Result_Scalar => Complex,
89 Left_Vector => Complex_Vector,
90 Result_Vector => Complex_Vector,
91 Operation => "*");
93 function "*" is new Vector_Scalar_Elementwise_Operation
94 (Left_Scalar => Complex,
95 Right_Scalar => Real'Base,
96 Result_Scalar => Complex,
97 Left_Vector => Complex_Vector,
98 Result_Vector => Complex_Vector,
99 Operation => "*");
101 function "*" is new Scalar_Vector_Elementwise_Operation
102 (Left_Scalar => Complex,
103 Right_Scalar => Complex,
104 Result_Scalar => Complex,
105 Right_Vector => Complex_Vector,
106 Result_Vector => Complex_Vector,
107 Operation => "*");
109 function "*" is new Scalar_Vector_Elementwise_Operation
110 (Left_Scalar => Real'Base,
111 Right_Scalar => Complex,
112 Result_Scalar => Complex,
113 Right_Vector => Complex_Vector,
114 Result_Vector => Complex_Vector,
115 Operation => "*");
117 function "*" is new Inner_Product
118 (Left_Scalar => Complex,
119 Right_Scalar => Real'Base,
120 Result_Scalar => Complex,
121 Left_Vector => Complex_Vector,
122 Right_Vector => Real_Vector,
123 Zero => (0.0, 0.0));
125 function "*" is new Inner_Product
126 (Left_Scalar => Real'Base,
127 Right_Scalar => Complex,
128 Result_Scalar => Complex,
129 Left_Vector => Real_Vector,
130 Right_Vector => Complex_Vector,
131 Zero => (0.0, 0.0));
133 function "*" is new Inner_Product
134 (Left_Scalar => Complex,
135 Right_Scalar => Complex,
136 Result_Scalar => Complex,
137 Left_Vector => Complex_Vector,
138 Right_Vector => Complex_Vector,
139 Zero => (0.0, 0.0));
141 function "*" is new Outer_Product
142 (Left_Scalar => Complex,
143 Right_Scalar => Complex,
144 Result_Scalar => Complex,
145 Left_Vector => Complex_Vector,
146 Right_Vector => Complex_Vector,
147 Matrix => Complex_Matrix);
149 function "*" is new Outer_Product
150 (Left_Scalar => Real'Base,
151 Right_Scalar => Complex,
152 Result_Scalar => Complex,
153 Left_Vector => Real_Vector,
154 Right_Vector => Complex_Vector,
155 Matrix => Complex_Matrix);
157 function "*" is new Outer_Product
158 (Left_Scalar => Complex,
159 Right_Scalar => Real'Base,
160 Result_Scalar => Complex,
161 Left_Vector => Complex_Vector,
162 Right_Vector => Real_Vector,
163 Matrix => Complex_Matrix);
165 function "*" is new Matrix_Scalar_Elementwise_Operation
166 (Left_Scalar => Complex,
167 Right_Scalar => Complex,
168 Result_Scalar => Complex,
169 Left_Matrix => Complex_Matrix,
170 Result_Matrix => Complex_Matrix,
171 Operation => "*");
173 function "*" is new Matrix_Scalar_Elementwise_Operation
174 (Left_Scalar => Complex,
175 Right_Scalar => Real'Base,
176 Result_Scalar => Complex,
177 Left_Matrix => Complex_Matrix,
178 Result_Matrix => Complex_Matrix,
179 Operation => "*");
181 function "*" is new Scalar_Matrix_Elementwise_Operation
182 (Left_Scalar => Complex,
183 Right_Scalar => Complex,
184 Result_Scalar => Complex,
185 Right_Matrix => Complex_Matrix,
186 Result_Matrix => Complex_Matrix,
187 Operation => "*");
189 function "*" is new Scalar_Matrix_Elementwise_Operation
190 (Left_Scalar => Real'Base,
191 Right_Scalar => Complex,
192 Result_Scalar => Complex,
193 Right_Matrix => Complex_Matrix,
194 Result_Matrix => Complex_Matrix,
195 Operation => "*");
197 function "*" is new Matrix_Vector_Product
198 (Left_Scalar => Real'Base,
199 Right_Scalar => Complex,
200 Result_Scalar => Complex,
201 Matrix => Real_Matrix,
202 Right_Vector => Complex_Vector,
203 Result_Vector => Complex_Vector,
204 Zero => (0.0, 0.0));
206 function "*" is new Matrix_Vector_Product
207 (Left_Scalar => Complex,
208 Right_Scalar => Real'Base,
209 Result_Scalar => Complex,
210 Matrix => Complex_Matrix,
211 Right_Vector => Real_Vector,
212 Result_Vector => Complex_Vector,
213 Zero => (0.0, 0.0));
215 function "*" is new Matrix_Vector_Product
216 (Left_Scalar => Complex,
217 Right_Scalar => Complex,
218 Result_Scalar => Complex,
219 Matrix => Complex_Matrix,
220 Right_Vector => Complex_Vector,
221 Result_Vector => Complex_Vector,
222 Zero => (0.0, 0.0));
224 function "*" is new Vector_Matrix_Product
225 (Left_Scalar => Real'Base,
226 Right_Scalar => Complex,
227 Result_Scalar => Complex,
228 Left_Vector => Real_Vector,
229 Matrix => Complex_Matrix,
230 Result_Vector => Complex_Vector,
231 Zero => (0.0, 0.0));
233 function "*" is new Vector_Matrix_Product
234 (Left_Scalar => Complex,
235 Right_Scalar => Real'Base,
236 Result_Scalar => Complex,
237 Left_Vector => Complex_Vector,
238 Matrix => Real_Matrix,
239 Result_Vector => Complex_Vector,
240 Zero => (0.0, 0.0));
242 function "*" is new Vector_Matrix_Product
243 (Left_Scalar => Complex,
244 Right_Scalar => Complex,
245 Result_Scalar => Complex,
246 Left_Vector => Complex_Vector,
247 Matrix => Complex_Matrix,
248 Result_Vector => Complex_Vector,
249 Zero => (0.0, 0.0));
251 function "*" is new Matrix_Matrix_Product
252 (Left_Scalar => Complex,
253 Right_Scalar => Complex,
254 Result_Scalar => Complex,
255 Left_Matrix => Complex_Matrix,
256 Right_Matrix => Complex_Matrix,
257 Result_Matrix => Complex_Matrix,
258 Zero => (0.0, 0.0));
260 function "*" is new Matrix_Matrix_Product
261 (Left_Scalar => Real'Base,
262 Right_Scalar => Complex,
263 Result_Scalar => Complex,
264 Left_Matrix => Real_Matrix,
265 Right_Matrix => Complex_Matrix,
266 Result_Matrix => Complex_Matrix,
267 Zero => (0.0, 0.0));
269 function "*" is new Matrix_Matrix_Product
270 (Left_Scalar => Complex,
271 Right_Scalar => Real'Base,
272 Result_Scalar => Complex,
273 Left_Matrix => Complex_Matrix,
274 Right_Matrix => Real_Matrix,
275 Result_Matrix => Complex_Matrix,
276 Zero => (0.0, 0.0));
278 ---------
279 -- "+" --
280 ---------
282 function "+" is new Vector_Elementwise_Operation
283 (X_Scalar => Complex,
284 Result_Scalar => Complex,
285 X_Vector => Complex_Vector,
286 Result_Vector => Complex_Vector,
287 Operation => "+");
289 function "+" is new Vector_Vector_Elementwise_Operation
290 (Left_Scalar => Complex,
291 Right_Scalar => Complex,
292 Result_Scalar => Complex,
293 Left_Vector => Complex_Vector,
294 Right_Vector => Complex_Vector,
295 Result_Vector => Complex_Vector,
296 Operation => "+");
298 function "+" is new Vector_Vector_Elementwise_Operation
299 (Left_Scalar => Real'Base,
300 Right_Scalar => Complex,
301 Result_Scalar => Complex,
302 Left_Vector => Real_Vector,
303 Right_Vector => Complex_Vector,
304 Result_Vector => Complex_Vector,
305 Operation => "+");
307 function "+" is new Vector_Vector_Elementwise_Operation
308 (Left_Scalar => Complex,
309 Right_Scalar => Real'Base,
310 Result_Scalar => Complex,
311 Left_Vector => Complex_Vector,
312 Right_Vector => Real_Vector,
313 Result_Vector => Complex_Vector,
314 Operation => "+");
316 function "+" is new Matrix_Elementwise_Operation
317 (X_Scalar => Complex,
318 Result_Scalar => Complex,
319 X_Matrix => Complex_Matrix,
320 Result_Matrix => Complex_Matrix,
321 Operation => "+");
323 function "+" is new Matrix_Matrix_Elementwise_Operation
324 (Left_Scalar => Complex,
325 Right_Scalar => Complex,
326 Result_Scalar => Complex,
327 Left_Matrix => Complex_Matrix,
328 Right_Matrix => Complex_Matrix,
329 Result_Matrix => Complex_Matrix,
330 Operation => "+");
332 function "+" is new Matrix_Matrix_Elementwise_Operation
333 (Left_Scalar => Real'Base,
334 Right_Scalar => Complex,
335 Result_Scalar => Complex,
336 Left_Matrix => Real_Matrix,
337 Right_Matrix => Complex_Matrix,
338 Result_Matrix => Complex_Matrix,
339 Operation => "+");
341 function "+" is new Matrix_Matrix_Elementwise_Operation
342 (Left_Scalar => Complex,
343 Right_Scalar => Real'Base,
344 Result_Scalar => Complex,
345 Left_Matrix => Complex_Matrix,
346 Right_Matrix => Real_Matrix,
347 Result_Matrix => Complex_Matrix,
348 Operation => "+");
350 ---------
351 -- "-" --
352 ---------
354 function "-" is new Vector_Elementwise_Operation
355 (X_Scalar => Complex,
356 Result_Scalar => Complex,
357 X_Vector => Complex_Vector,
358 Result_Vector => Complex_Vector,
359 Operation => "-");
361 function "-" is new Vector_Vector_Elementwise_Operation
362 (Left_Scalar => Complex,
363 Right_Scalar => Complex,
364 Result_Scalar => Complex,
365 Left_Vector => Complex_Vector,
366 Right_Vector => Complex_Vector,
367 Result_Vector => Complex_Vector,
368 Operation => "-");
370 function "-" is new Vector_Vector_Elementwise_Operation
371 (Left_Scalar => Real'Base,
372 Right_Scalar => Complex,
373 Result_Scalar => Complex,
374 Left_Vector => Real_Vector,
375 Right_Vector => Complex_Vector,
376 Result_Vector => Complex_Vector,
377 Operation => "-");
379 function "-" is new Vector_Vector_Elementwise_Operation
380 (Left_Scalar => Complex,
381 Right_Scalar => Real'Base,
382 Result_Scalar => Complex,
383 Left_Vector => Complex_Vector,
384 Right_Vector => Real_Vector,
385 Result_Vector => Complex_Vector,
386 Operation => "-");
388 function "-" is new Matrix_Elementwise_Operation
389 (X_Scalar => Complex,
390 Result_Scalar => Complex,
391 X_Matrix => Complex_Matrix,
392 Result_Matrix => Complex_Matrix,
393 Operation => "-");
395 function "-" is new Matrix_Matrix_Elementwise_Operation
396 (Left_Scalar => Complex,
397 Right_Scalar => Complex,
398 Result_Scalar => Complex,
399 Left_Matrix => Complex_Matrix,
400 Right_Matrix => Complex_Matrix,
401 Result_Matrix => Complex_Matrix,
402 Operation => "-");
404 function "-" is new Matrix_Matrix_Elementwise_Operation
405 (Left_Scalar => Real'Base,
406 Right_Scalar => Complex,
407 Result_Scalar => Complex,
408 Left_Matrix => Real_Matrix,
409 Right_Matrix => Complex_Matrix,
410 Result_Matrix => Complex_Matrix,
411 Operation => "-");
413 function "-" is new Matrix_Matrix_Elementwise_Operation
414 (Left_Scalar => Complex,
415 Right_Scalar => Real'Base,
416 Result_Scalar => Complex,
417 Left_Matrix => Complex_Matrix,
418 Right_Matrix => Real_Matrix,
419 Result_Matrix => Complex_Matrix,
420 Operation => "-");
422 ---------
423 -- "/" --
424 ---------
426 function "/" is new Vector_Scalar_Elementwise_Operation
427 (Left_Scalar => Complex,
428 Right_Scalar => Complex,
429 Result_Scalar => Complex,
430 Left_Vector => Complex_Vector,
431 Result_Vector => Complex_Vector,
432 Operation => "/");
434 function "/" is new Vector_Scalar_Elementwise_Operation
435 (Left_Scalar => Complex,
436 Right_Scalar => Real'Base,
437 Result_Scalar => Complex,
438 Left_Vector => Complex_Vector,
439 Result_Vector => Complex_Vector,
440 Operation => "/");
442 function "/" is new Matrix_Scalar_Elementwise_Operation
443 (Left_Scalar => Complex,
444 Right_Scalar => Complex,
445 Result_Scalar => Complex,
446 Left_Matrix => Complex_Matrix,
447 Result_Matrix => Complex_Matrix,
448 Operation => "/");
450 function "/" is new Matrix_Scalar_Elementwise_Operation
451 (Left_Scalar => Complex,
452 Right_Scalar => Real'Base,
453 Result_Scalar => Complex,
454 Left_Matrix => Complex_Matrix,
455 Result_Matrix => Complex_Matrix,
456 Operation => "/");
458 -----------
459 -- "abs" --
460 -----------
462 function "abs" is new L2_Norm
463 (X_Scalar => Complex,
464 Result_Real => Real'Base,
465 X_Vector => Complex_Vector);
467 --------------
468 -- Argument --
469 --------------
471 function Argument is new Vector_Elementwise_Operation
472 (X_Scalar => Complex,
473 Result_Scalar => Real'Base,
474 X_Vector => Complex_Vector,
475 Result_Vector => Real_Vector,
476 Operation => Argument);
478 function Argument is new Vector_Scalar_Elementwise_Operation
479 (Left_Scalar => Complex,
480 Right_Scalar => Real'Base,
481 Result_Scalar => Real'Base,
482 Left_Vector => Complex_Vector,
483 Result_Vector => Real_Vector,
484 Operation => Argument);
486 function Argument is new Matrix_Elementwise_Operation
487 (X_Scalar => Complex,
488 Result_Scalar => Real'Base,
489 X_Matrix => Complex_Matrix,
490 Result_Matrix => Real_Matrix,
491 Operation => Argument);
493 function Argument is new Matrix_Scalar_Elementwise_Operation
494 (Left_Scalar => Complex,
495 Right_Scalar => Real'Base,
496 Result_Scalar => Real'Base,
497 Left_Matrix => Complex_Matrix,
498 Result_Matrix => Real_Matrix,
499 Operation => Argument);
501 ----------------------------
502 -- Compose_From_Cartesian --
503 ----------------------------
505 function Compose_From_Cartesian is new Vector_Elementwise_Operation
506 (X_Scalar => Real'Base,
507 Result_Scalar => Complex,
508 X_Vector => Real_Vector,
509 Result_Vector => Complex_Vector,
510 Operation => Compose_From_Cartesian);
512 function Compose_From_Cartesian is
513 new Vector_Vector_Elementwise_Operation
514 (Left_Scalar => Real'Base,
515 Right_Scalar => Real'Base,
516 Result_Scalar => Complex,
517 Left_Vector => Real_Vector,
518 Right_Vector => Real_Vector,
519 Result_Vector => Complex_Vector,
520 Operation => Compose_From_Cartesian);
522 function Compose_From_Cartesian is new Matrix_Elementwise_Operation
523 (X_Scalar => Real'Base,
524 Result_Scalar => Complex,
525 X_Matrix => Real_Matrix,
526 Result_Matrix => Complex_Matrix,
527 Operation => Compose_From_Cartesian);
529 function Compose_From_Cartesian is
530 new Matrix_Matrix_Elementwise_Operation
531 (Left_Scalar => Real'Base,
532 Right_Scalar => Real'Base,
533 Result_Scalar => Complex,
534 Left_Matrix => Real_Matrix,
535 Right_Matrix => Real_Matrix,
536 Result_Matrix => Complex_Matrix,
537 Operation => Compose_From_Cartesian);
539 ------------------------
540 -- Compose_From_Polar --
541 ------------------------
543 function Compose_From_Polar is
544 new Vector_Vector_Elementwise_Operation
545 (Left_Scalar => Real'Base,
546 Right_Scalar => Real'Base,
547 Result_Scalar => Complex,
548 Left_Vector => Real_Vector,
549 Right_Vector => Real_Vector,
550 Result_Vector => Complex_Vector,
551 Operation => Compose_From_Polar);
553 function Compose_From_Polar is
554 new Vector_Vector_Scalar_Elementwise_Operation
555 (X_Scalar => Real'Base,
556 Y_Scalar => Real'Base,
557 Z_Scalar => Real'Base,
558 Result_Scalar => Complex,
559 X_Vector => Real_Vector,
560 Y_Vector => Real_Vector,
561 Result_Vector => Complex_Vector,
562 Operation => Compose_From_Polar);
564 function Compose_From_Polar is
565 new Matrix_Matrix_Elementwise_Operation
566 (Left_Scalar => Real'Base,
567 Right_Scalar => Real'Base,
568 Result_Scalar => Complex,
569 Left_Matrix => Real_Matrix,
570 Right_Matrix => Real_Matrix,
571 Result_Matrix => Complex_Matrix,
572 Operation => Compose_From_Polar);
574 function Compose_From_Polar is
575 new Matrix_Matrix_Scalar_Elementwise_Operation
576 (X_Scalar => Real'Base,
577 Y_Scalar => Real'Base,
578 Z_Scalar => Real'Base,
579 Result_Scalar => Complex,
580 X_Matrix => Real_Matrix,
581 Y_Matrix => Real_Matrix,
582 Result_Matrix => Complex_Matrix,
583 Operation => Compose_From_Polar);
585 ---------------
586 -- Conjugate --
587 ---------------
589 function Conjugate is new Vector_Elementwise_Operation
590 (X_Scalar => Complex,
591 Result_Scalar => Complex,
592 X_Vector => Complex_Vector,
593 Result_Vector => Complex_Vector,
594 Operation => Conjugate);
596 function Conjugate is new Matrix_Elementwise_Operation
597 (X_Scalar => Complex,
598 Result_Scalar => Complex,
599 X_Matrix => Complex_Matrix,
600 Result_Matrix => Complex_Matrix,
601 Operation => Conjugate);
603 --------
604 -- Im --
605 --------
607 function Im is new Vector_Elementwise_Operation
608 (X_Scalar => Complex,
609 Result_Scalar => Real'Base,
610 X_Vector => Complex_Vector,
611 Result_Vector => Real_Vector,
612 Operation => Im);
614 function Im is new Matrix_Elementwise_Operation
615 (X_Scalar => Complex,
616 Result_Scalar => Real'Base,
617 X_Matrix => Complex_Matrix,
618 Result_Matrix => Real_Matrix,
619 Operation => Im);
621 -------------
622 -- Modulus --
623 -------------
625 function Modulus is new Vector_Elementwise_Operation
626 (X_Scalar => Complex,
627 Result_Scalar => Real'Base,
628 X_Vector => Complex_Vector,
629 Result_Vector => Real_Vector,
630 Operation => Modulus);
632 function Modulus is new Matrix_Elementwise_Operation
633 (X_Scalar => Complex,
634 Result_Scalar => Real'Base,
635 X_Matrix => Complex_Matrix,
636 Result_Matrix => Real_Matrix,
637 Operation => Modulus);
639 --------
640 -- Re --
641 --------
643 function Re is new Vector_Elementwise_Operation
644 (X_Scalar => Complex,
645 Result_Scalar => Real'Base,
646 X_Vector => Complex_Vector,
647 Result_Vector => Real_Vector,
648 Operation => Re);
650 function Re is new Matrix_Elementwise_Operation
651 (X_Scalar => Complex,
652 Result_Scalar => Real'Base,
653 X_Matrix => Complex_Matrix,
654 Result_Matrix => Real_Matrix,
655 Operation => Re);
657 ------------
658 -- Set_Im --
659 ------------
661 procedure Set_Im is new Update_Vector_With_Vector
662 (X_Scalar => Complex,
663 Y_Scalar => Real'Base,
664 X_Vector => Complex_Vector,
665 Y_Vector => Real_Vector,
666 Update => Set_Im);
668 procedure Set_Im is new Update_Matrix_With_Matrix
669 (X_Scalar => Complex,
670 Y_Scalar => Real'Base,
671 X_Matrix => Complex_Matrix,
672 Y_Matrix => Real_Matrix,
673 Update => Set_Im);
675 ------------
676 -- Set_Re --
677 ------------
679 procedure Set_Re is new Update_Vector_With_Vector
680 (X_Scalar => Complex,
681 Y_Scalar => Real'Base,
682 X_Vector => Complex_Vector,
683 Y_Vector => Real_Vector,
684 Update => Set_Re);
686 procedure Set_Re is new Update_Matrix_With_Matrix
687 (X_Scalar => Complex,
688 Y_Scalar => Real'Base,
689 X_Matrix => Complex_Matrix,
690 Y_Matrix => Real_Matrix,
691 Update => Set_Re);
693 -----------
694 -- Solve --
695 -----------
697 function Solve is
698 new Matrix_Vector_Solution (Complex, Complex_Vector, Complex_Matrix);
700 function Solve is
701 new Matrix_Matrix_Solution (Complex, Complex_Matrix);
703 -----------------
704 -- Unit_Matrix --
705 -----------------
707 function Unit_Matrix is new System.Generic_Array_Operations.Unit_Matrix
708 (Scalar => Complex,
709 Matrix => Complex_Matrix,
710 Zero => (0.0, 0.0),
711 One => (1.0, 0.0));
713 function Unit_Vector is new System.Generic_Array_Operations.Unit_Vector
714 (Scalar => Complex,
715 Vector => Complex_Vector,
716 Zero => (0.0, 0.0),
717 One => (1.0, 0.0));
718 end Instantiations;
720 ---------
721 -- "*" --
722 ---------
724 function "*"
725 (Left : Complex_Vector;
726 Right : Complex_Vector) return Complex
727 renames Instantiations."*";
729 function "*"
730 (Left : Real_Vector;
731 Right : Complex_Vector) return Complex
732 renames Instantiations."*";
734 function "*"
735 (Left : Complex_Vector;
736 Right : Real_Vector) return Complex
737 renames Instantiations."*";
739 function "*"
740 (Left : Complex;
741 Right : Complex_Vector) return Complex_Vector
742 renames Instantiations."*";
744 function "*"
745 (Left : Complex_Vector;
746 Right : Complex) return Complex_Vector
747 renames Instantiations."*";
749 function "*"
750 (Left : Real'Base;
751 Right : Complex_Vector) return Complex_Vector
752 renames Instantiations."*";
754 function "*"
755 (Left : Complex_Vector;
756 Right : Real'Base) return Complex_Vector
757 renames Instantiations."*";
759 function "*"
760 (Left : Complex_Matrix;
761 Right : Complex_Matrix) return Complex_Matrix
762 renames Instantiations."*";
764 function "*"
765 (Left : Complex_Vector;
766 Right : Complex_Vector) return Complex_Matrix
767 renames Instantiations."*";
769 function "*"
770 (Left : Complex_Vector;
771 Right : Complex_Matrix) return Complex_Vector
772 renames Instantiations."*";
774 function "*"
775 (Left : Complex_Matrix;
776 Right : Complex_Vector) return Complex_Vector
777 renames Instantiations."*";
779 function "*"
780 (Left : Real_Matrix;
781 Right : Complex_Matrix) return Complex_Matrix
782 renames Instantiations."*";
784 function "*"
785 (Left : Complex_Matrix;
786 Right : Real_Matrix) return Complex_Matrix
787 renames Instantiations."*";
789 function "*"
790 (Left : Real_Vector;
791 Right : Complex_Vector) return Complex_Matrix
792 renames Instantiations."*";
794 function "*"
795 (Left : Complex_Vector;
796 Right : Real_Vector) return Complex_Matrix
797 renames Instantiations."*";
799 function "*"
800 (Left : Real_Vector;
801 Right : Complex_Matrix) return Complex_Vector
802 renames Instantiations."*";
804 function "*"
805 (Left : Complex_Vector;
806 Right : Real_Matrix) return Complex_Vector
807 renames Instantiations."*";
809 function "*"
810 (Left : Real_Matrix;
811 Right : Complex_Vector) return Complex_Vector
812 renames Instantiations."*";
814 function "*"
815 (Left : Complex_Matrix;
816 Right : Real_Vector) return Complex_Vector
817 renames Instantiations."*";
819 function "*"
820 (Left : Complex;
821 Right : Complex_Matrix) return Complex_Matrix
822 renames Instantiations."*";
824 function "*"
825 (Left : Complex_Matrix;
826 Right : Complex) return Complex_Matrix
827 renames Instantiations."*";
829 function "*"
830 (Left : Real'Base;
831 Right : Complex_Matrix) return Complex_Matrix
832 renames Instantiations."*";
834 function "*"
835 (Left : Complex_Matrix;
836 Right : Real'Base) return Complex_Matrix
837 renames Instantiations."*";
839 ---------
840 -- "+" --
841 ---------
843 function "+" (Right : Complex_Vector) return Complex_Vector
844 renames Instantiations."+";
846 function "+"
847 (Left : Complex_Vector;
848 Right : Complex_Vector) return Complex_Vector
849 renames Instantiations."+";
851 function "+"
852 (Left : Real_Vector;
853 Right : Complex_Vector) return Complex_Vector
854 renames Instantiations."+";
856 function "+"
857 (Left : Complex_Vector;
858 Right : Real_Vector) return Complex_Vector
859 renames Instantiations."+";
861 function "+" (Right : Complex_Matrix) return Complex_Matrix
862 renames Instantiations."+";
864 function "+"
865 (Left : Complex_Matrix;
866 Right : Complex_Matrix) return Complex_Matrix
867 renames Instantiations."+";
869 function "+"
870 (Left : Real_Matrix;
871 Right : Complex_Matrix) return Complex_Matrix
872 renames Instantiations."+";
874 function "+"
875 (Left : Complex_Matrix;
876 Right : Real_Matrix) return Complex_Matrix
877 renames Instantiations."+";
879 ---------
880 -- "-" --
881 ---------
883 function "-"
884 (Right : Complex_Vector) return Complex_Vector
885 renames Instantiations."-";
887 function "-"
888 (Left : Complex_Vector;
889 Right : Complex_Vector) return Complex_Vector
890 renames Instantiations."-";
892 function "-"
893 (Left : Real_Vector;
894 Right : Complex_Vector) return Complex_Vector
895 renames Instantiations."-";
897 function "-"
898 (Left : Complex_Vector;
899 Right : Real_Vector) return Complex_Vector
900 renames Instantiations."-";
902 function "-" (Right : Complex_Matrix) return Complex_Matrix
903 renames Instantiations."-";
905 function "-"
906 (Left : Complex_Matrix;
907 Right : Complex_Matrix) return Complex_Matrix
908 renames Instantiations."-";
910 function "-"
911 (Left : Real_Matrix;
912 Right : Complex_Matrix) return Complex_Matrix
913 renames Instantiations."-";
915 function "-"
916 (Left : Complex_Matrix;
917 Right : Real_Matrix) return Complex_Matrix
918 renames Instantiations."-";
920 ---------
921 -- "/" --
922 ---------
924 function "/"
925 (Left : Complex_Vector;
926 Right : Complex) return Complex_Vector
927 renames Instantiations."/";
929 function "/"
930 (Left : Complex_Vector;
931 Right : Real'Base) return Complex_Vector
932 renames Instantiations."/";
934 function "/"
935 (Left : Complex_Matrix;
936 Right : Complex) return Complex_Matrix
937 renames Instantiations."/";
939 function "/"
940 (Left : Complex_Matrix;
941 Right : Real'Base) return Complex_Matrix
942 renames Instantiations."/";
944 -----------
945 -- "abs" --
946 -----------
948 function "abs" (Right : Complex_Vector) return Real'Base
949 renames Instantiations."abs";
951 --------------
952 -- Argument --
953 --------------
955 function Argument (X : Complex_Vector) return Real_Vector
956 renames Instantiations.Argument;
958 function Argument
959 (X : Complex_Vector;
960 Cycle : Real'Base) return Real_Vector
961 renames Instantiations.Argument;
963 function Argument (X : Complex_Matrix) return Real_Matrix
964 renames Instantiations.Argument;
966 function Argument
967 (X : Complex_Matrix;
968 Cycle : Real'Base) return Real_Matrix
969 renames Instantiations.Argument;
971 ----------------------------
972 -- Compose_From_Cartesian --
973 ----------------------------
975 function Compose_From_Cartesian (Re : Real_Vector) return Complex_Vector
976 renames Instantiations.Compose_From_Cartesian;
978 function Compose_From_Cartesian
979 (Re : Real_Vector;
980 Im : Real_Vector) return Complex_Vector
981 renames Instantiations.Compose_From_Cartesian;
983 function Compose_From_Cartesian (Re : Real_Matrix) return Complex_Matrix
984 renames Instantiations.Compose_From_Cartesian;
986 function Compose_From_Cartesian
987 (Re : Real_Matrix;
988 Im : Real_Matrix) return Complex_Matrix
989 renames Instantiations.Compose_From_Cartesian;
991 ------------------------
992 -- Compose_From_Polar --
993 ------------------------
995 function Compose_From_Polar
996 (Modulus : Real_Vector;
997 Argument : Real_Vector) return Complex_Vector
998 renames Instantiations.Compose_From_Polar;
1000 function Compose_From_Polar
1001 (Modulus : Real_Vector;
1002 Argument : Real_Vector;
1003 Cycle : Real'Base) return Complex_Vector
1004 renames Instantiations.Compose_From_Polar;
1006 function Compose_From_Polar
1007 (Modulus : Real_Matrix;
1008 Argument : Real_Matrix) return Complex_Matrix
1009 renames Instantiations.Compose_From_Polar;
1011 function Compose_From_Polar
1012 (Modulus : Real_Matrix;
1013 Argument : Real_Matrix;
1014 Cycle : Real'Base) return Complex_Matrix
1015 renames Instantiations.Compose_From_Polar;
1017 ---------------
1018 -- Conjugate --
1019 ---------------
1021 function Conjugate (X : Complex_Vector) return Complex_Vector
1022 renames Instantiations.Conjugate;
1024 function Conjugate (X : Complex_Matrix) return Complex_Matrix
1025 renames Instantiations.Conjugate;
1027 -----------------
1028 -- Determinant --
1029 -----------------
1031 function Determinant (A : Complex_Matrix) return Complex is
1032 M : Complex_Matrix := A;
1033 B : Complex_Matrix (A'Range (1), 1 .. 0);
1034 R : Complex;
1035 begin
1036 Forward_Eliminate (M, B, R);
1037 return R;
1038 end Determinant;
1040 -----------------
1041 -- Eigensystem --
1042 -----------------
1044 procedure Eigensystem
1045 (A : Complex_Matrix;
1046 Values : out Real_Vector;
1047 Vectors : out Complex_Matrix)
1049 N : constant Natural := Length (A);
1051 -- For a Hermitian matrix C, we convert the eigenvalue problem to a
1052 -- real symmetric one: if C = A + i * B, then the (N, N) complex
1053 -- eigenvalue problem:
1054 -- (A + i * B) * (u + i * v) = Lambda * (u + i * v)
1056 -- is equivalent to the (2 * N, 2 * N) real eigenvalue problem:
1057 -- [ A, B ] [ u ] = Lambda * [ u ]
1058 -- [ -B, A ] [ v ] [ v ]
1060 -- Note that the (2 * N, 2 * N) matrix above is symmetric, as
1061 -- Transpose (A) = A and Transpose (B) = -B if C is Hermitian.
1063 -- We solve this eigensystem using the real-valued algorithms. The final
1064 -- result will have every eigenvalue twice, so in the sorted output we
1065 -- just pick every second value, with associated eigenvector u + i * v.
1067 M : Real_Matrix (1 .. 2 * N, 1 .. 2 * N);
1068 Vals : Real_Vector (1 .. 2 * N);
1069 Vecs : Real_Matrix (1 .. 2 * N, 1 .. 2 * N);
1071 begin
1072 for J in 1 .. N loop
1073 for K in 1 .. N loop
1074 declare
1075 C : constant Complex :=
1076 (A (A'First (1) + (J - 1), A'First (2) + (K - 1)));
1077 begin
1078 M (J, K) := Re (C);
1079 M (J + N, K + N) := Re (C);
1080 M (J + N, K) := Im (C);
1081 M (J, K + N) := -Im (C);
1082 end;
1083 end loop;
1084 end loop;
1086 Eigensystem (M, Vals, Vecs);
1088 for J in 1 .. N loop
1089 declare
1090 Col : constant Integer := Values'First + (J - 1);
1091 begin
1092 Values (Col) := Vals (2 * J);
1094 for K in 1 .. N loop
1095 declare
1096 Row : constant Integer := Vectors'First (2) + (K - 1);
1097 begin
1098 Vectors (Row, Col)
1099 := (Vecs (J * 2, Col), Vecs (J * 2, Col + N));
1100 end;
1101 end loop;
1102 end;
1103 end loop;
1104 end Eigensystem;
1106 -----------------
1107 -- Eigenvalues --
1108 -----------------
1110 function Eigenvalues (A : Complex_Matrix) return Real_Vector is
1111 -- See Eigensystem for a description of the algorithm
1113 N : constant Natural := Length (A);
1114 R : Real_Vector (A'Range (1));
1116 M : Real_Matrix (1 .. 2 * N, 1 .. 2 * N);
1117 Vals : Real_Vector (1 .. 2 * N);
1118 begin
1119 for J in 1 .. N loop
1120 for K in 1 .. N loop
1121 declare
1122 C : constant Complex :=
1123 (A (A'First (1) + (J - 1), A'First (2) + (K - 1)));
1124 begin
1125 M (J, K) := Re (C);
1126 M (J + N, K + N) := Re (C);
1127 M (J + N, K) := Im (C);
1128 M (J, K + N) := -Im (C);
1129 end;
1130 end loop;
1131 end loop;
1133 Vals := Eigenvalues (M);
1135 for J in 1 .. N loop
1136 R (A'First (1) + (J - 1)) := Vals (2 * J);
1137 end loop;
1139 return R;
1140 end Eigenvalues;
1142 --------
1143 -- Im --
1144 --------
1146 function Im (X : Complex_Vector) return Real_Vector
1147 renames Instantiations.Im;
1149 function Im (X : Complex_Matrix) return Real_Matrix
1150 renames Instantiations.Im;
1152 -------------
1153 -- Inverse --
1154 -------------
1156 function Inverse (A : Complex_Matrix) return Complex_Matrix is
1157 (Solve (A, Unit_Matrix (Length (A))));
1159 -------------
1160 -- Modulus --
1161 -------------
1163 function Modulus (X : Complex_Vector) return Real_Vector
1164 renames Instantiations.Modulus;
1166 function Modulus (X : Complex_Matrix) return Real_Matrix
1167 renames Instantiations.Modulus;
1169 --------
1170 -- Re --
1171 --------
1173 function Re (X : Complex_Vector) return Real_Vector
1174 renames Instantiations.Re;
1176 function Re (X : Complex_Matrix) return Real_Matrix
1177 renames Instantiations.Re;
1179 ------------
1180 -- Set_Im --
1181 ------------
1183 procedure Set_Im
1184 (X : in out Complex_Matrix;
1185 Im : Real_Matrix)
1186 renames Instantiations.Set_Im;
1188 procedure Set_Im
1189 (X : in out Complex_Vector;
1190 Im : Real_Vector)
1191 renames Instantiations.Set_Im;
1193 ------------
1194 -- Set_Re --
1195 ------------
1197 procedure Set_Re
1198 (X : in out Complex_Matrix;
1199 Re : Real_Matrix)
1200 renames Instantiations.Set_Re;
1202 procedure Set_Re
1203 (X : in out Complex_Vector;
1204 Re : Real_Vector)
1205 renames Instantiations.Set_Re;
1207 -----------
1208 -- Solve --
1209 -----------
1211 function Solve
1212 (A : Complex_Matrix;
1213 X : Complex_Vector) return Complex_Vector
1214 renames Instantiations.Solve;
1216 function Solve
1217 (A : Complex_Matrix;
1218 X : Complex_Matrix) return Complex_Matrix
1219 renames Instantiations.Solve;
1221 ---------------
1222 -- Transpose --
1223 ---------------
1225 function Transpose
1226 (X : Complex_Matrix) return Complex_Matrix
1228 R : Complex_Matrix (X'Range (2), X'Range (1));
1229 begin
1230 Transpose (X, R);
1231 return R;
1232 end Transpose;
1234 -----------------
1235 -- Unit_Matrix --
1236 -----------------
1238 function Unit_Matrix
1239 (Order : Positive;
1240 First_1 : Integer := 1;
1241 First_2 : Integer := 1) return Complex_Matrix
1242 renames Instantiations.Unit_Matrix;
1244 -----------------
1245 -- Unit_Vector --
1246 -----------------
1248 function Unit_Vector
1249 (Index : Integer;
1250 Order : Positive;
1251 First : Integer := 1) return Complex_Vector
1252 renames Instantiations.Unit_Vector;
1254 end Ada.Numerics.Generic_Complex_Arrays;