PR rtl-optimization/82913
[official-gcc.git] / gcc / ada / libgnat / a-ngcoar.adb
blobcf01dcdef806889adf57e41ca35c35b289016b9b
1 ------------------------------------------------------------------------------
2 -- --
3 -- GNAT COMPILER COMPONENTS --
4 -- --
5 -- ADA.NUMERICS.GENERIC_COMPLEX_ARRAYS --
6 -- --
7 -- B o d y --
8 -- --
9 -- Copyright (C) 2006-2017, 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;
34 package body Ada.Numerics.Generic_Complex_Arrays is
36 -- Operations that are defined in terms of operations on the type Real,
37 -- such as addition, subtraction and scaling, are computed in the canonical
38 -- way looping over all elements.
40 package Ops renames System.Generic_Array_Operations;
42 subtype Real is Real_Arrays.Real;
43 -- Work around visibility bug ???
45 function Is_Non_Zero (X : Complex) return Boolean is (X /= (0.0, 0.0));
46 -- Needed by Back_Substitute
48 procedure Back_Substitute is new Ops.Back_Substitute
49 (Scalar => Complex,
50 Matrix => Complex_Matrix,
51 Is_Non_Zero => Is_Non_Zero);
53 procedure Forward_Eliminate is new Ops.Forward_Eliminate
54 (Scalar => Complex,
55 Real => Real'Base,
56 Matrix => Complex_Matrix,
57 Zero => (0.0, 0.0),
58 One => (1.0, 0.0));
60 procedure Transpose is new Ops.Transpose
61 (Scalar => Complex,
62 Matrix => Complex_Matrix);
64 -- Helper function that raises a Constraint_Error is the argument is
65 -- not a square matrix, and otherwise returns its length.
67 function Length is new Square_Matrix_Length (Complex, Complex_Matrix);
69 -- Instant a generic square root implementation here, in order to avoid
70 -- instantiating a complete copy of Generic_Elementary_Functions.
71 -- Speed of the square root is not a big concern here.
73 function Sqrt is new Ops.Sqrt (Real'Base);
75 -- Instantiating the following subprograms directly would lead to
76 -- name clashes, so use a local package.
78 package Instantiations is
80 ---------
81 -- "*" --
82 ---------
84 function "*" is new Vector_Scalar_Elementwise_Operation
85 (Left_Scalar => Complex,
86 Right_Scalar => Complex,
87 Result_Scalar => Complex,
88 Left_Vector => Complex_Vector,
89 Result_Vector => Complex_Vector,
90 Operation => "*");
92 function "*" is new Vector_Scalar_Elementwise_Operation
93 (Left_Scalar => Complex,
94 Right_Scalar => Real'Base,
95 Result_Scalar => Complex,
96 Left_Vector => Complex_Vector,
97 Result_Vector => Complex_Vector,
98 Operation => "*");
100 function "*" is new Scalar_Vector_Elementwise_Operation
101 (Left_Scalar => Complex,
102 Right_Scalar => Complex,
103 Result_Scalar => Complex,
104 Right_Vector => Complex_Vector,
105 Result_Vector => Complex_Vector,
106 Operation => "*");
108 function "*" is new Scalar_Vector_Elementwise_Operation
109 (Left_Scalar => Real'Base,
110 Right_Scalar => Complex,
111 Result_Scalar => Complex,
112 Right_Vector => Complex_Vector,
113 Result_Vector => Complex_Vector,
114 Operation => "*");
116 function "*" is new Inner_Product
117 (Left_Scalar => Complex,
118 Right_Scalar => Real'Base,
119 Result_Scalar => Complex,
120 Left_Vector => Complex_Vector,
121 Right_Vector => Real_Vector,
122 Zero => (0.0, 0.0));
124 function "*" is new Inner_Product
125 (Left_Scalar => Real'Base,
126 Right_Scalar => Complex,
127 Result_Scalar => Complex,
128 Left_Vector => Real_Vector,
129 Right_Vector => Complex_Vector,
130 Zero => (0.0, 0.0));
132 function "*" is new Inner_Product
133 (Left_Scalar => Complex,
134 Right_Scalar => Complex,
135 Result_Scalar => Complex,
136 Left_Vector => Complex_Vector,
137 Right_Vector => Complex_Vector,
138 Zero => (0.0, 0.0));
140 function "*" is new Outer_Product
141 (Left_Scalar => Complex,
142 Right_Scalar => Complex,
143 Result_Scalar => Complex,
144 Left_Vector => Complex_Vector,
145 Right_Vector => Complex_Vector,
146 Matrix => Complex_Matrix);
148 function "*" is new Outer_Product
149 (Left_Scalar => Real'Base,
150 Right_Scalar => Complex,
151 Result_Scalar => Complex,
152 Left_Vector => Real_Vector,
153 Right_Vector => Complex_Vector,
154 Matrix => Complex_Matrix);
156 function "*" is new Outer_Product
157 (Left_Scalar => Complex,
158 Right_Scalar => Real'Base,
159 Result_Scalar => Complex,
160 Left_Vector => Complex_Vector,
161 Right_Vector => Real_Vector,
162 Matrix => Complex_Matrix);
164 function "*" is new Matrix_Scalar_Elementwise_Operation
165 (Left_Scalar => Complex,
166 Right_Scalar => Complex,
167 Result_Scalar => Complex,
168 Left_Matrix => Complex_Matrix,
169 Result_Matrix => Complex_Matrix,
170 Operation => "*");
172 function "*" is new Matrix_Scalar_Elementwise_Operation
173 (Left_Scalar => Complex,
174 Right_Scalar => Real'Base,
175 Result_Scalar => Complex,
176 Left_Matrix => Complex_Matrix,
177 Result_Matrix => Complex_Matrix,
178 Operation => "*");
180 function "*" is new Scalar_Matrix_Elementwise_Operation
181 (Left_Scalar => Complex,
182 Right_Scalar => Complex,
183 Result_Scalar => Complex,
184 Right_Matrix => Complex_Matrix,
185 Result_Matrix => Complex_Matrix,
186 Operation => "*");
188 function "*" is new Scalar_Matrix_Elementwise_Operation
189 (Left_Scalar => Real'Base,
190 Right_Scalar => Complex,
191 Result_Scalar => Complex,
192 Right_Matrix => Complex_Matrix,
193 Result_Matrix => Complex_Matrix,
194 Operation => "*");
196 function "*" is new Matrix_Vector_Product
197 (Left_Scalar => Real'Base,
198 Right_Scalar => Complex,
199 Result_Scalar => Complex,
200 Matrix => Real_Matrix,
201 Right_Vector => Complex_Vector,
202 Result_Vector => Complex_Vector,
203 Zero => (0.0, 0.0));
205 function "*" is new Matrix_Vector_Product
206 (Left_Scalar => Complex,
207 Right_Scalar => Real'Base,
208 Result_Scalar => Complex,
209 Matrix => Complex_Matrix,
210 Right_Vector => Real_Vector,
211 Result_Vector => Complex_Vector,
212 Zero => (0.0, 0.0));
214 function "*" is new Matrix_Vector_Product
215 (Left_Scalar => Complex,
216 Right_Scalar => Complex,
217 Result_Scalar => Complex,
218 Matrix => Complex_Matrix,
219 Right_Vector => Complex_Vector,
220 Result_Vector => Complex_Vector,
221 Zero => (0.0, 0.0));
223 function "*" is new Vector_Matrix_Product
224 (Left_Scalar => Real'Base,
225 Right_Scalar => Complex,
226 Result_Scalar => Complex,
227 Left_Vector => Real_Vector,
228 Matrix => Complex_Matrix,
229 Result_Vector => Complex_Vector,
230 Zero => (0.0, 0.0));
232 function "*" is new Vector_Matrix_Product
233 (Left_Scalar => Complex,
234 Right_Scalar => Real'Base,
235 Result_Scalar => Complex,
236 Left_Vector => Complex_Vector,
237 Matrix => Real_Matrix,
238 Result_Vector => Complex_Vector,
239 Zero => (0.0, 0.0));
241 function "*" is new Vector_Matrix_Product
242 (Left_Scalar => Complex,
243 Right_Scalar => Complex,
244 Result_Scalar => Complex,
245 Left_Vector => Complex_Vector,
246 Matrix => Complex_Matrix,
247 Result_Vector => Complex_Vector,
248 Zero => (0.0, 0.0));
250 function "*" is new Matrix_Matrix_Product
251 (Left_Scalar => Complex,
252 Right_Scalar => Complex,
253 Result_Scalar => Complex,
254 Left_Matrix => Complex_Matrix,
255 Right_Matrix => Complex_Matrix,
256 Result_Matrix => Complex_Matrix,
257 Zero => (0.0, 0.0));
259 function "*" is new Matrix_Matrix_Product
260 (Left_Scalar => Real'Base,
261 Right_Scalar => Complex,
262 Result_Scalar => Complex,
263 Left_Matrix => Real_Matrix,
264 Right_Matrix => Complex_Matrix,
265 Result_Matrix => Complex_Matrix,
266 Zero => (0.0, 0.0));
268 function "*" is new Matrix_Matrix_Product
269 (Left_Scalar => Complex,
270 Right_Scalar => Real'Base,
271 Result_Scalar => Complex,
272 Left_Matrix => Complex_Matrix,
273 Right_Matrix => Real_Matrix,
274 Result_Matrix => Complex_Matrix,
275 Zero => (0.0, 0.0));
277 ---------
278 -- "+" --
279 ---------
281 function "+" is new Vector_Elementwise_Operation
282 (X_Scalar => Complex,
283 Result_Scalar => Complex,
284 X_Vector => Complex_Vector,
285 Result_Vector => Complex_Vector,
286 Operation => "+");
288 function "+" is new Vector_Vector_Elementwise_Operation
289 (Left_Scalar => Complex,
290 Right_Scalar => Complex,
291 Result_Scalar => Complex,
292 Left_Vector => Complex_Vector,
293 Right_Vector => Complex_Vector,
294 Result_Vector => Complex_Vector,
295 Operation => "+");
297 function "+" is new Vector_Vector_Elementwise_Operation
298 (Left_Scalar => Real'Base,
299 Right_Scalar => Complex,
300 Result_Scalar => Complex,
301 Left_Vector => Real_Vector,
302 Right_Vector => Complex_Vector,
303 Result_Vector => Complex_Vector,
304 Operation => "+");
306 function "+" is new Vector_Vector_Elementwise_Operation
307 (Left_Scalar => Complex,
308 Right_Scalar => Real'Base,
309 Result_Scalar => Complex,
310 Left_Vector => Complex_Vector,
311 Right_Vector => Real_Vector,
312 Result_Vector => Complex_Vector,
313 Operation => "+");
315 function "+" is new Matrix_Elementwise_Operation
316 (X_Scalar => Complex,
317 Result_Scalar => Complex,
318 X_Matrix => Complex_Matrix,
319 Result_Matrix => Complex_Matrix,
320 Operation => "+");
322 function "+" is new Matrix_Matrix_Elementwise_Operation
323 (Left_Scalar => Complex,
324 Right_Scalar => Complex,
325 Result_Scalar => Complex,
326 Left_Matrix => Complex_Matrix,
327 Right_Matrix => Complex_Matrix,
328 Result_Matrix => Complex_Matrix,
329 Operation => "+");
331 function "+" is new Matrix_Matrix_Elementwise_Operation
332 (Left_Scalar => Real'Base,
333 Right_Scalar => Complex,
334 Result_Scalar => Complex,
335 Left_Matrix => Real_Matrix,
336 Right_Matrix => Complex_Matrix,
337 Result_Matrix => Complex_Matrix,
338 Operation => "+");
340 function "+" is new Matrix_Matrix_Elementwise_Operation
341 (Left_Scalar => Complex,
342 Right_Scalar => Real'Base,
343 Result_Scalar => Complex,
344 Left_Matrix => Complex_Matrix,
345 Right_Matrix => Real_Matrix,
346 Result_Matrix => Complex_Matrix,
347 Operation => "+");
349 ---------
350 -- "-" --
351 ---------
353 function "-" is new Vector_Elementwise_Operation
354 (X_Scalar => Complex,
355 Result_Scalar => Complex,
356 X_Vector => Complex_Vector,
357 Result_Vector => Complex_Vector,
358 Operation => "-");
360 function "-" is new Vector_Vector_Elementwise_Operation
361 (Left_Scalar => Complex,
362 Right_Scalar => Complex,
363 Result_Scalar => Complex,
364 Left_Vector => Complex_Vector,
365 Right_Vector => Complex_Vector,
366 Result_Vector => Complex_Vector,
367 Operation => "-");
369 function "-" is new Vector_Vector_Elementwise_Operation
370 (Left_Scalar => Real'Base,
371 Right_Scalar => Complex,
372 Result_Scalar => Complex,
373 Left_Vector => Real_Vector,
374 Right_Vector => Complex_Vector,
375 Result_Vector => Complex_Vector,
376 Operation => "-");
378 function "-" is new Vector_Vector_Elementwise_Operation
379 (Left_Scalar => Complex,
380 Right_Scalar => Real'Base,
381 Result_Scalar => Complex,
382 Left_Vector => Complex_Vector,
383 Right_Vector => Real_Vector,
384 Result_Vector => Complex_Vector,
385 Operation => "-");
387 function "-" is new Matrix_Elementwise_Operation
388 (X_Scalar => Complex,
389 Result_Scalar => Complex,
390 X_Matrix => Complex_Matrix,
391 Result_Matrix => Complex_Matrix,
392 Operation => "-");
394 function "-" is new Matrix_Matrix_Elementwise_Operation
395 (Left_Scalar => Complex,
396 Right_Scalar => Complex,
397 Result_Scalar => Complex,
398 Left_Matrix => Complex_Matrix,
399 Right_Matrix => Complex_Matrix,
400 Result_Matrix => Complex_Matrix,
401 Operation => "-");
403 function "-" is new Matrix_Matrix_Elementwise_Operation
404 (Left_Scalar => Real'Base,
405 Right_Scalar => Complex,
406 Result_Scalar => Complex,
407 Left_Matrix => Real_Matrix,
408 Right_Matrix => Complex_Matrix,
409 Result_Matrix => Complex_Matrix,
410 Operation => "-");
412 function "-" is new Matrix_Matrix_Elementwise_Operation
413 (Left_Scalar => Complex,
414 Right_Scalar => Real'Base,
415 Result_Scalar => Complex,
416 Left_Matrix => Complex_Matrix,
417 Right_Matrix => Real_Matrix,
418 Result_Matrix => Complex_Matrix,
419 Operation => "-");
421 ---------
422 -- "/" --
423 ---------
425 function "/" is new Vector_Scalar_Elementwise_Operation
426 (Left_Scalar => Complex,
427 Right_Scalar => Complex,
428 Result_Scalar => Complex,
429 Left_Vector => Complex_Vector,
430 Result_Vector => Complex_Vector,
431 Operation => "/");
433 function "/" is new Vector_Scalar_Elementwise_Operation
434 (Left_Scalar => Complex,
435 Right_Scalar => Real'Base,
436 Result_Scalar => Complex,
437 Left_Vector => Complex_Vector,
438 Result_Vector => Complex_Vector,
439 Operation => "/");
441 function "/" is new Matrix_Scalar_Elementwise_Operation
442 (Left_Scalar => Complex,
443 Right_Scalar => Complex,
444 Result_Scalar => Complex,
445 Left_Matrix => Complex_Matrix,
446 Result_Matrix => Complex_Matrix,
447 Operation => "/");
449 function "/" is new Matrix_Scalar_Elementwise_Operation
450 (Left_Scalar => Complex,
451 Right_Scalar => Real'Base,
452 Result_Scalar => Complex,
453 Left_Matrix => Complex_Matrix,
454 Result_Matrix => Complex_Matrix,
455 Operation => "/");
457 -----------
458 -- "abs" --
459 -----------
461 function "abs" is new L2_Norm
462 (X_Scalar => Complex,
463 Result_Real => Real'Base,
464 X_Vector => Complex_Vector);
466 --------------
467 -- Argument --
468 --------------
470 function Argument is new Vector_Elementwise_Operation
471 (X_Scalar => Complex,
472 Result_Scalar => Real'Base,
473 X_Vector => Complex_Vector,
474 Result_Vector => Real_Vector,
475 Operation => Argument);
477 function Argument is new Vector_Scalar_Elementwise_Operation
478 (Left_Scalar => Complex,
479 Right_Scalar => Real'Base,
480 Result_Scalar => Real'Base,
481 Left_Vector => Complex_Vector,
482 Result_Vector => Real_Vector,
483 Operation => Argument);
485 function Argument is new Matrix_Elementwise_Operation
486 (X_Scalar => Complex,
487 Result_Scalar => Real'Base,
488 X_Matrix => Complex_Matrix,
489 Result_Matrix => Real_Matrix,
490 Operation => Argument);
492 function Argument is new Matrix_Scalar_Elementwise_Operation
493 (Left_Scalar => Complex,
494 Right_Scalar => Real'Base,
495 Result_Scalar => Real'Base,
496 Left_Matrix => Complex_Matrix,
497 Result_Matrix => Real_Matrix,
498 Operation => Argument);
500 ----------------------------
501 -- Compose_From_Cartesian --
502 ----------------------------
504 function Compose_From_Cartesian is new Vector_Elementwise_Operation
505 (X_Scalar => Real'Base,
506 Result_Scalar => Complex,
507 X_Vector => Real_Vector,
508 Result_Vector => Complex_Vector,
509 Operation => Compose_From_Cartesian);
511 function Compose_From_Cartesian is
512 new Vector_Vector_Elementwise_Operation
513 (Left_Scalar => Real'Base,
514 Right_Scalar => Real'Base,
515 Result_Scalar => Complex,
516 Left_Vector => Real_Vector,
517 Right_Vector => Real_Vector,
518 Result_Vector => Complex_Vector,
519 Operation => Compose_From_Cartesian);
521 function Compose_From_Cartesian is new Matrix_Elementwise_Operation
522 (X_Scalar => Real'Base,
523 Result_Scalar => Complex,
524 X_Matrix => Real_Matrix,
525 Result_Matrix => Complex_Matrix,
526 Operation => Compose_From_Cartesian);
528 function Compose_From_Cartesian is
529 new Matrix_Matrix_Elementwise_Operation
530 (Left_Scalar => Real'Base,
531 Right_Scalar => Real'Base,
532 Result_Scalar => Complex,
533 Left_Matrix => Real_Matrix,
534 Right_Matrix => Real_Matrix,
535 Result_Matrix => Complex_Matrix,
536 Operation => Compose_From_Cartesian);
538 ------------------------
539 -- Compose_From_Polar --
540 ------------------------
542 function Compose_From_Polar is
543 new Vector_Vector_Elementwise_Operation
544 (Left_Scalar => Real'Base,
545 Right_Scalar => Real'Base,
546 Result_Scalar => Complex,
547 Left_Vector => Real_Vector,
548 Right_Vector => Real_Vector,
549 Result_Vector => Complex_Vector,
550 Operation => Compose_From_Polar);
552 function Compose_From_Polar is
553 new Vector_Vector_Scalar_Elementwise_Operation
554 (X_Scalar => Real'Base,
555 Y_Scalar => Real'Base,
556 Z_Scalar => Real'Base,
557 Result_Scalar => Complex,
558 X_Vector => Real_Vector,
559 Y_Vector => Real_Vector,
560 Result_Vector => Complex_Vector,
561 Operation => Compose_From_Polar);
563 function Compose_From_Polar is
564 new Matrix_Matrix_Elementwise_Operation
565 (Left_Scalar => Real'Base,
566 Right_Scalar => Real'Base,
567 Result_Scalar => Complex,
568 Left_Matrix => Real_Matrix,
569 Right_Matrix => Real_Matrix,
570 Result_Matrix => Complex_Matrix,
571 Operation => Compose_From_Polar);
573 function Compose_From_Polar is
574 new Matrix_Matrix_Scalar_Elementwise_Operation
575 (X_Scalar => Real'Base,
576 Y_Scalar => Real'Base,
577 Z_Scalar => Real'Base,
578 Result_Scalar => Complex,
579 X_Matrix => Real_Matrix,
580 Y_Matrix => Real_Matrix,
581 Result_Matrix => Complex_Matrix,
582 Operation => Compose_From_Polar);
584 ---------------
585 -- Conjugate --
586 ---------------
588 function Conjugate is new Vector_Elementwise_Operation
589 (X_Scalar => Complex,
590 Result_Scalar => Complex,
591 X_Vector => Complex_Vector,
592 Result_Vector => Complex_Vector,
593 Operation => Conjugate);
595 function Conjugate is new Matrix_Elementwise_Operation
596 (X_Scalar => Complex,
597 Result_Scalar => Complex,
598 X_Matrix => Complex_Matrix,
599 Result_Matrix => Complex_Matrix,
600 Operation => Conjugate);
602 --------
603 -- Im --
604 --------
606 function Im is new Vector_Elementwise_Operation
607 (X_Scalar => Complex,
608 Result_Scalar => Real'Base,
609 X_Vector => Complex_Vector,
610 Result_Vector => Real_Vector,
611 Operation => Im);
613 function Im is new Matrix_Elementwise_Operation
614 (X_Scalar => Complex,
615 Result_Scalar => Real'Base,
616 X_Matrix => Complex_Matrix,
617 Result_Matrix => Real_Matrix,
618 Operation => Im);
620 -------------
621 -- Modulus --
622 -------------
624 function Modulus is new Vector_Elementwise_Operation
625 (X_Scalar => Complex,
626 Result_Scalar => Real'Base,
627 X_Vector => Complex_Vector,
628 Result_Vector => Real_Vector,
629 Operation => Modulus);
631 function Modulus is new Matrix_Elementwise_Operation
632 (X_Scalar => Complex,
633 Result_Scalar => Real'Base,
634 X_Matrix => Complex_Matrix,
635 Result_Matrix => Real_Matrix,
636 Operation => Modulus);
638 --------
639 -- Re --
640 --------
642 function Re is new Vector_Elementwise_Operation
643 (X_Scalar => Complex,
644 Result_Scalar => Real'Base,
645 X_Vector => Complex_Vector,
646 Result_Vector => Real_Vector,
647 Operation => Re);
649 function Re is new Matrix_Elementwise_Operation
650 (X_Scalar => Complex,
651 Result_Scalar => Real'Base,
652 X_Matrix => Complex_Matrix,
653 Result_Matrix => Real_Matrix,
654 Operation => Re);
656 ------------
657 -- Set_Im --
658 ------------
660 procedure Set_Im is new Update_Vector_With_Vector
661 (X_Scalar => Complex,
662 Y_Scalar => Real'Base,
663 X_Vector => Complex_Vector,
664 Y_Vector => Real_Vector,
665 Update => Set_Im);
667 procedure Set_Im is new Update_Matrix_With_Matrix
668 (X_Scalar => Complex,
669 Y_Scalar => Real'Base,
670 X_Matrix => Complex_Matrix,
671 Y_Matrix => Real_Matrix,
672 Update => Set_Im);
674 ------------
675 -- Set_Re --
676 ------------
678 procedure Set_Re is new Update_Vector_With_Vector
679 (X_Scalar => Complex,
680 Y_Scalar => Real'Base,
681 X_Vector => Complex_Vector,
682 Y_Vector => Real_Vector,
683 Update => Set_Re);
685 procedure Set_Re is new Update_Matrix_With_Matrix
686 (X_Scalar => Complex,
687 Y_Scalar => Real'Base,
688 X_Matrix => Complex_Matrix,
689 Y_Matrix => Real_Matrix,
690 Update => Set_Re);
692 -----------
693 -- Solve --
694 -----------
696 function Solve is new Matrix_Vector_Solution
697 (Complex, (0.0, 0.0), Complex_Vector, Complex_Matrix);
699 function Solve is new Matrix_Matrix_Solution
700 (Complex, (0.0, 0.0), Complex_Matrix);
702 -----------------
703 -- Unit_Matrix --
704 -----------------
706 function Unit_Matrix is new System.Generic_Array_Operations.Unit_Matrix
707 (Scalar => Complex,
708 Matrix => Complex_Matrix,
709 Zero => (0.0, 0.0),
710 One => (1.0, 0.0));
712 function Unit_Vector is new System.Generic_Array_Operations.Unit_Vector
713 (Scalar => Complex,
714 Vector => Complex_Vector,
715 Zero => (0.0, 0.0),
716 One => (1.0, 0.0));
717 end Instantiations;
719 ---------
720 -- "*" --
721 ---------
723 function "*"
724 (Left : Complex_Vector;
725 Right : Complex_Vector) return Complex
726 renames Instantiations."*";
728 function "*"
729 (Left : Real_Vector;
730 Right : Complex_Vector) return Complex
731 renames Instantiations."*";
733 function "*"
734 (Left : Complex_Vector;
735 Right : Real_Vector) return Complex
736 renames Instantiations."*";
738 function "*"
739 (Left : Complex;
740 Right : Complex_Vector) return Complex_Vector
741 renames Instantiations."*";
743 function "*"
744 (Left : Complex_Vector;
745 Right : Complex) return Complex_Vector
746 renames Instantiations."*";
748 function "*"
749 (Left : Real'Base;
750 Right : Complex_Vector) return Complex_Vector
751 renames Instantiations."*";
753 function "*"
754 (Left : Complex_Vector;
755 Right : Real'Base) return Complex_Vector
756 renames Instantiations."*";
758 function "*"
759 (Left : Complex_Matrix;
760 Right : Complex_Matrix) return Complex_Matrix
761 renames Instantiations."*";
763 function "*"
764 (Left : Complex_Vector;
765 Right : Complex_Vector) return Complex_Matrix
766 renames Instantiations."*";
768 function "*"
769 (Left : Complex_Vector;
770 Right : Complex_Matrix) return Complex_Vector
771 renames Instantiations."*";
773 function "*"
774 (Left : Complex_Matrix;
775 Right : Complex_Vector) return Complex_Vector
776 renames Instantiations."*";
778 function "*"
779 (Left : Real_Matrix;
780 Right : Complex_Matrix) return Complex_Matrix
781 renames Instantiations."*";
783 function "*"
784 (Left : Complex_Matrix;
785 Right : Real_Matrix) return Complex_Matrix
786 renames Instantiations."*";
788 function "*"
789 (Left : Real_Vector;
790 Right : Complex_Vector) return Complex_Matrix
791 renames Instantiations."*";
793 function "*"
794 (Left : Complex_Vector;
795 Right : Real_Vector) return Complex_Matrix
796 renames Instantiations."*";
798 function "*"
799 (Left : Real_Vector;
800 Right : Complex_Matrix) return Complex_Vector
801 renames Instantiations."*";
803 function "*"
804 (Left : Complex_Vector;
805 Right : Real_Matrix) return Complex_Vector
806 renames Instantiations."*";
808 function "*"
809 (Left : Real_Matrix;
810 Right : Complex_Vector) return Complex_Vector
811 renames Instantiations."*";
813 function "*"
814 (Left : Complex_Matrix;
815 Right : Real_Vector) return Complex_Vector
816 renames Instantiations."*";
818 function "*"
819 (Left : Complex;
820 Right : Complex_Matrix) return Complex_Matrix
821 renames Instantiations."*";
823 function "*"
824 (Left : Complex_Matrix;
825 Right : Complex) return Complex_Matrix
826 renames Instantiations."*";
828 function "*"
829 (Left : Real'Base;
830 Right : Complex_Matrix) return Complex_Matrix
831 renames Instantiations."*";
833 function "*"
834 (Left : Complex_Matrix;
835 Right : Real'Base) return Complex_Matrix
836 renames Instantiations."*";
838 ---------
839 -- "+" --
840 ---------
842 function "+" (Right : Complex_Vector) return Complex_Vector
843 renames Instantiations."+";
845 function "+"
846 (Left : Complex_Vector;
847 Right : Complex_Vector) return Complex_Vector
848 renames Instantiations."+";
850 function "+"
851 (Left : Real_Vector;
852 Right : Complex_Vector) return Complex_Vector
853 renames Instantiations."+";
855 function "+"
856 (Left : Complex_Vector;
857 Right : Real_Vector) return Complex_Vector
858 renames Instantiations."+";
860 function "+" (Right : Complex_Matrix) return Complex_Matrix
861 renames Instantiations."+";
863 function "+"
864 (Left : Complex_Matrix;
865 Right : Complex_Matrix) return Complex_Matrix
866 renames Instantiations."+";
868 function "+"
869 (Left : Real_Matrix;
870 Right : Complex_Matrix) return Complex_Matrix
871 renames Instantiations."+";
873 function "+"
874 (Left : Complex_Matrix;
875 Right : Real_Matrix) return Complex_Matrix
876 renames Instantiations."+";
878 ---------
879 -- "-" --
880 ---------
882 function "-"
883 (Right : Complex_Vector) return Complex_Vector
884 renames Instantiations."-";
886 function "-"
887 (Left : Complex_Vector;
888 Right : Complex_Vector) return Complex_Vector
889 renames Instantiations."-";
891 function "-"
892 (Left : Real_Vector;
893 Right : Complex_Vector) return Complex_Vector
894 renames Instantiations."-";
896 function "-"
897 (Left : Complex_Vector;
898 Right : Real_Vector) return Complex_Vector
899 renames Instantiations."-";
901 function "-" (Right : Complex_Matrix) return Complex_Matrix
902 renames Instantiations."-";
904 function "-"
905 (Left : Complex_Matrix;
906 Right : Complex_Matrix) return Complex_Matrix
907 renames Instantiations."-";
909 function "-"
910 (Left : Real_Matrix;
911 Right : Complex_Matrix) return Complex_Matrix
912 renames Instantiations."-";
914 function "-"
915 (Left : Complex_Matrix;
916 Right : Real_Matrix) return Complex_Matrix
917 renames Instantiations."-";
919 ---------
920 -- "/" --
921 ---------
923 function "/"
924 (Left : Complex_Vector;
925 Right : Complex) return Complex_Vector
926 renames Instantiations."/";
928 function "/"
929 (Left : Complex_Vector;
930 Right : Real'Base) return Complex_Vector
931 renames Instantiations."/";
933 function "/"
934 (Left : Complex_Matrix;
935 Right : Complex) return Complex_Matrix
936 renames Instantiations."/";
938 function "/"
939 (Left : Complex_Matrix;
940 Right : Real'Base) return Complex_Matrix
941 renames Instantiations."/";
943 -----------
944 -- "abs" --
945 -----------
947 function "abs" (Right : Complex_Vector) return Real'Base
948 renames Instantiations."abs";
950 --------------
951 -- Argument --
952 --------------
954 function Argument (X : Complex_Vector) return Real_Vector
955 renames Instantiations.Argument;
957 function Argument
958 (X : Complex_Vector;
959 Cycle : Real'Base) return Real_Vector
960 renames Instantiations.Argument;
962 function Argument (X : Complex_Matrix) return Real_Matrix
963 renames Instantiations.Argument;
965 function Argument
966 (X : Complex_Matrix;
967 Cycle : Real'Base) return Real_Matrix
968 renames Instantiations.Argument;
970 ----------------------------
971 -- Compose_From_Cartesian --
972 ----------------------------
974 function Compose_From_Cartesian (Re : Real_Vector) return Complex_Vector
975 renames Instantiations.Compose_From_Cartesian;
977 function Compose_From_Cartesian
978 (Re : Real_Vector;
979 Im : Real_Vector) return Complex_Vector
980 renames Instantiations.Compose_From_Cartesian;
982 function Compose_From_Cartesian (Re : Real_Matrix) return Complex_Matrix
983 renames Instantiations.Compose_From_Cartesian;
985 function Compose_From_Cartesian
986 (Re : Real_Matrix;
987 Im : Real_Matrix) return Complex_Matrix
988 renames Instantiations.Compose_From_Cartesian;
990 ------------------------
991 -- Compose_From_Polar --
992 ------------------------
994 function Compose_From_Polar
995 (Modulus : Real_Vector;
996 Argument : Real_Vector) return Complex_Vector
997 renames Instantiations.Compose_From_Polar;
999 function Compose_From_Polar
1000 (Modulus : Real_Vector;
1001 Argument : Real_Vector;
1002 Cycle : Real'Base) return Complex_Vector
1003 renames Instantiations.Compose_From_Polar;
1005 function Compose_From_Polar
1006 (Modulus : Real_Matrix;
1007 Argument : Real_Matrix) return Complex_Matrix
1008 renames Instantiations.Compose_From_Polar;
1010 function Compose_From_Polar
1011 (Modulus : Real_Matrix;
1012 Argument : Real_Matrix;
1013 Cycle : Real'Base) return Complex_Matrix
1014 renames Instantiations.Compose_From_Polar;
1016 ---------------
1017 -- Conjugate --
1018 ---------------
1020 function Conjugate (X : Complex_Vector) return Complex_Vector
1021 renames Instantiations.Conjugate;
1023 function Conjugate (X : Complex_Matrix) return Complex_Matrix
1024 renames Instantiations.Conjugate;
1026 -----------------
1027 -- Determinant --
1028 -----------------
1030 function Determinant (A : Complex_Matrix) return Complex is
1031 M : Complex_Matrix := A;
1032 B : Complex_Matrix (A'Range (1), 1 .. 0);
1033 R : Complex;
1034 begin
1035 Forward_Eliminate (M, B, R);
1036 return R;
1037 end Determinant;
1039 -----------------
1040 -- Eigensystem --
1041 -----------------
1043 procedure Eigensystem
1044 (A : Complex_Matrix;
1045 Values : out Real_Vector;
1046 Vectors : out Complex_Matrix)
1048 N : constant Natural := Length (A);
1050 -- For a Hermitian matrix C, we convert the eigenvalue problem to a
1051 -- real symmetric one: if C = A + i * B, then the (N, N) complex
1052 -- eigenvalue problem:
1053 -- (A + i * B) * (u + i * v) = Lambda * (u + i * v)
1055 -- is equivalent to the (2 * N, 2 * N) real eigenvalue problem:
1056 -- [ A, B ] [ u ] = Lambda * [ u ]
1057 -- [ -B, A ] [ v ] [ v ]
1059 -- Note that the (2 * N, 2 * N) matrix above is symmetric, as
1060 -- Transpose (A) = A and Transpose (B) = -B if C is Hermitian.
1062 -- We solve this eigensystem using the real-valued algorithms. The final
1063 -- result will have every eigenvalue twice, so in the sorted output we
1064 -- just pick every second value, with associated eigenvector u + i * v.
1066 M : Real_Matrix (1 .. 2 * N, 1 .. 2 * N);
1067 Vals : Real_Vector (1 .. 2 * N);
1068 Vecs : Real_Matrix (1 .. 2 * N, 1 .. 2 * N);
1070 begin
1071 for J in 1 .. N loop
1072 for K in 1 .. N loop
1073 declare
1074 C : constant Complex :=
1075 (A (A'First (1) + (J - 1), A'First (2) + (K - 1)));
1076 begin
1077 M (J, K) := Re (C);
1078 M (J + N, K + N) := Re (C);
1079 M (J + N, K) := Im (C);
1080 M (J, K + N) := -Im (C);
1081 end;
1082 end loop;
1083 end loop;
1085 Eigensystem (M, Vals, Vecs);
1087 for J in 1 .. N loop
1088 declare
1089 Col : constant Integer := Values'First + (J - 1);
1090 begin
1091 Values (Col) := Vals (2 * J);
1093 for K in 1 .. N loop
1094 declare
1095 Row : constant Integer := Vectors'First (2) + (K - 1);
1096 begin
1097 Vectors (Row, Col)
1098 := (Vecs (J * 2, Col), Vecs (J * 2, Col + N));
1099 end;
1100 end loop;
1101 end;
1102 end loop;
1103 end Eigensystem;
1105 -----------------
1106 -- Eigenvalues --
1107 -----------------
1109 function Eigenvalues (A : Complex_Matrix) return Real_Vector is
1110 -- See Eigensystem for a description of the algorithm
1112 N : constant Natural := Length (A);
1113 R : Real_Vector (A'Range (1));
1115 M : Real_Matrix (1 .. 2 * N, 1 .. 2 * N);
1116 Vals : Real_Vector (1 .. 2 * N);
1117 begin
1118 for J in 1 .. N loop
1119 for K in 1 .. N loop
1120 declare
1121 C : constant Complex :=
1122 (A (A'First (1) + (J - 1), A'First (2) + (K - 1)));
1123 begin
1124 M (J, K) := Re (C);
1125 M (J + N, K + N) := Re (C);
1126 M (J + N, K) := Im (C);
1127 M (J, K + N) := -Im (C);
1128 end;
1129 end loop;
1130 end loop;
1132 Vals := Eigenvalues (M);
1134 for J in 1 .. N loop
1135 R (A'First (1) + (J - 1)) := Vals (2 * J);
1136 end loop;
1138 return R;
1139 end Eigenvalues;
1141 --------
1142 -- Im --
1143 --------
1145 function Im (X : Complex_Vector) return Real_Vector
1146 renames Instantiations.Im;
1148 function Im (X : Complex_Matrix) return Real_Matrix
1149 renames Instantiations.Im;
1151 -------------
1152 -- Inverse --
1153 -------------
1155 function Inverse (A : Complex_Matrix) return Complex_Matrix is
1156 (Solve (A, Unit_Matrix (Length (A),
1157 First_1 => A'First (2),
1158 First_2 => A'First (1))));
1160 -------------
1161 -- Modulus --
1162 -------------
1164 function Modulus (X : Complex_Vector) return Real_Vector
1165 renames Instantiations.Modulus;
1167 function Modulus (X : Complex_Matrix) return Real_Matrix
1168 renames Instantiations.Modulus;
1170 --------
1171 -- Re --
1172 --------
1174 function Re (X : Complex_Vector) return Real_Vector
1175 renames Instantiations.Re;
1177 function Re (X : Complex_Matrix) return Real_Matrix
1178 renames Instantiations.Re;
1180 ------------
1181 -- Set_Im --
1182 ------------
1184 procedure Set_Im
1185 (X : in out Complex_Matrix;
1186 Im : Real_Matrix)
1187 renames Instantiations.Set_Im;
1189 procedure Set_Im
1190 (X : in out Complex_Vector;
1191 Im : Real_Vector)
1192 renames Instantiations.Set_Im;
1194 ------------
1195 -- Set_Re --
1196 ------------
1198 procedure Set_Re
1199 (X : in out Complex_Matrix;
1200 Re : Real_Matrix)
1201 renames Instantiations.Set_Re;
1203 procedure Set_Re
1204 (X : in out Complex_Vector;
1205 Re : Real_Vector)
1206 renames Instantiations.Set_Re;
1208 -----------
1209 -- Solve --
1210 -----------
1212 function Solve
1213 (A : Complex_Matrix;
1214 X : Complex_Vector) return Complex_Vector
1215 renames Instantiations.Solve;
1217 function Solve
1218 (A : Complex_Matrix;
1219 X : Complex_Matrix) return Complex_Matrix
1220 renames Instantiations.Solve;
1222 ---------------
1223 -- Transpose --
1224 ---------------
1226 function Transpose
1227 (X : Complex_Matrix) return Complex_Matrix
1229 R : Complex_Matrix (X'Range (2), X'Range (1));
1230 begin
1231 Transpose (X, R);
1232 return R;
1233 end Transpose;
1235 -----------------
1236 -- Unit_Matrix --
1237 -----------------
1239 function Unit_Matrix
1240 (Order : Positive;
1241 First_1 : Integer := 1;
1242 First_2 : Integer := 1) return Complex_Matrix
1243 renames Instantiations.Unit_Matrix;
1245 -----------------
1246 -- Unit_Vector --
1247 -----------------
1249 function Unit_Vector
1250 (Index : Integer;
1251 Order : Positive;
1252 First : Integer := 1) return Complex_Vector
1253 renames Instantiations.Unit_Vector;
1255 end Ada.Numerics.Generic_Complex_Arrays;