File : a-ngcoar.adb
1 ------------------------------------------------------------------------------
2 -- --
3 -- GNAT COMPILER COMPONENTS --
4 -- --
5 -- ADA.NUMERICS.GENERIC_COMPLEX_ARRAYS --
6 -- --
7 -- B o d y --
8 -- --
9 -- Copyright (C) 2006-2016, 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 -- --
19 -- --
20 -- --
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 ------------------------------------------------------------------------------
31
32 with System.Generic_Array_Operations; use System.Generic_Array_Operations;
33
34 package body Ada.Numerics.Generic_Complex_Arrays is
35
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.
39
40 package Ops renames System.Generic_Array_Operations;
41
42 subtype Real is Real_Arrays.Real;
43 -- Work around visibility bug ???
44
45 function Is_Non_Zero (X : Complex) return Boolean is (X /= (0.0, 0.0));
46 -- Needed by Back_Substitute
47
48 procedure Back_Substitute is new Ops.Back_Substitute
49 (Scalar => Complex,
50 Matrix => Complex_Matrix,
51 Is_Non_Zero => Is_Non_Zero);
52
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));
59
60 procedure Transpose is new Ops.Transpose
61 (Scalar => Complex,
62 Matrix => Complex_Matrix);
63
64 -- Helper function that raises a Constraint_Error is the argument is
65 -- not a square matrix, and otherwise returns its length.
66
67 function Length is new Square_Matrix_Length (Complex, Complex_Matrix);
68
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.
72
73 function Sqrt is new Ops.Sqrt (Real'Base);
74
75 -- Instantiating the following subprograms directly would lead to
76 -- name clashes, so use a local package.
77
78 package Instantiations is
79
80 ---------
81 -- "*" --
82 ---------
83
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 => "*");
91
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 => "*");
99
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 => "*");
107
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 => "*");
115
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));
123
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));
131
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));
139
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);
147
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);
155
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);
163
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 => "*");
171
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 => "*");
179
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 => "*");
187
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 => "*");
195
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));
204
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));
213
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));
222
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));
231
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));
240
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));
249
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));
258
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));
267
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));
276
277 ---------
278 -- "+" --
279 ---------
280
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 => "+");
287
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 => "+");
296
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 => "+");
305
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 => "+");
314
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 => "+");
321
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 => "+");
330
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 => "+");
339
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 => "+");
348
349 ---------
350 -- "-" --
351 ---------
352
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 => "-");
359
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 => "-");
368
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 => "-");
377
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 => "-");
386
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 => "-");
393
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 => "-");
402
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 => "-");
411
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 => "-");
420
421 ---------
422 -- "/" --
423 ---------
424
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 => "/");
432
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 => "/");
440
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 => "/");
448
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 => "/");
456
457 -----------
458 -- "abs" --
459 -----------
460
461 function "abs" is new L2_Norm
462 (X_Scalar => Complex,
463 Result_Real => Real'Base,
464 X_Vector => Complex_Vector);
465
466 --------------
467 -- Argument --
468 --------------
469
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);
476
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);
484
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);
491
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);
499
500 ----------------------------
501 -- Compose_From_Cartesian --
502 ----------------------------
503
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);
510
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);
520
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);
527
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);
537
538 ------------------------
539 -- Compose_From_Polar --
540 ------------------------
541
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);
551
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);
562
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);
572
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);
583
584 ---------------
585 -- Conjugate --
586 ---------------
587
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);
594
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);
601
602 --------
603 -- Im --
604 --------
605
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);
612
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);
619
620 -------------
621 -- Modulus --
622 -------------
623
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);
630
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);
637
638 --------
639 -- Re --
640 --------
641
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);
648
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);
655
656 ------------
657 -- Set_Im --
658 ------------
659
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);
666
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);
673
674 ------------
675 -- Set_Re --
676 ------------
677
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);
684
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);
691
692 -----------
693 -- Solve --
694 -----------
695
696 function Solve is new Matrix_Vector_Solution
697 (Complex, (0.0, 0.0), Complex_Vector, Complex_Matrix);
698
699 function Solve is new Matrix_Matrix_Solution
700 (Complex, (0.0, 0.0), Complex_Matrix);
701
702 -----------------
703 -- Unit_Matrix --
704 -----------------
705
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));
711
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;
718
719 ---------
720 -- "*" --
721 ---------
722
723 function "*"
724 (Left : Complex_Vector;
725 Right : Complex_Vector) return Complex
726 renames Instantiations."*";
727
728 function "*"
729 (Left : Real_Vector;
730 Right : Complex_Vector) return Complex
731 renames Instantiations."*";
732
733 function "*"
734 (Left : Complex_Vector;
735 Right : Real_Vector) return Complex
736 renames Instantiations."*";
737
738 function "*"
739 (Left : Complex;
740 Right : Complex_Vector) return Complex_Vector
741 renames Instantiations."*";
742
743 function "*"
744 (Left : Complex_Vector;
745 Right : Complex) return Complex_Vector
746 renames Instantiations."*";
747
748 function "*"
749 (Left : Real'Base;
750 Right : Complex_Vector) return Complex_Vector
751 renames Instantiations."*";
752
753 function "*"
754 (Left : Complex_Vector;
755 Right : Real'Base) return Complex_Vector
756 renames Instantiations."*";
757
758 function "*"
759 (Left : Complex_Matrix;
760 Right : Complex_Matrix) return Complex_Matrix
761 renames Instantiations."*";
762
763 function "*"
764 (Left : Complex_Vector;
765 Right : Complex_Vector) return Complex_Matrix
766 renames Instantiations."*";
767
768 function "*"
769 (Left : Complex_Vector;
770 Right : Complex_Matrix) return Complex_Vector
771 renames Instantiations."*";
772
773 function "*"
774 (Left : Complex_Matrix;
775 Right : Complex_Vector) return Complex_Vector
776 renames Instantiations."*";
777
778 function "*"
779 (Left : Real_Matrix;
780 Right : Complex_Matrix) return Complex_Matrix
781 renames Instantiations."*";
782
783 function "*"
784 (Left : Complex_Matrix;
785 Right : Real_Matrix) return Complex_Matrix
786 renames Instantiations."*";
787
788 function "*"
789 (Left : Real_Vector;
790 Right : Complex_Vector) return Complex_Matrix
791 renames Instantiations."*";
792
793 function "*"
794 (Left : Complex_Vector;
795 Right : Real_Vector) return Complex_Matrix
796 renames Instantiations."*";
797
798 function "*"
799 (Left : Real_Vector;
800 Right : Complex_Matrix) return Complex_Vector
801 renames Instantiations."*";
802
803 function "*"
804 (Left : Complex_Vector;
805 Right : Real_Matrix) return Complex_Vector
806 renames Instantiations."*";
807
808 function "*"
809 (Left : Real_Matrix;
810 Right : Complex_Vector) return Complex_Vector
811 renames Instantiations."*";
812
813 function "*"
814 (Left : Complex_Matrix;
815 Right : Real_Vector) return Complex_Vector
816 renames Instantiations."*";
817
818 function "*"
819 (Left : Complex;
820 Right : Complex_Matrix) return Complex_Matrix
821 renames Instantiations."*";
822
823 function "*"
824 (Left : Complex_Matrix;
825 Right : Complex) return Complex_Matrix
826 renames Instantiations."*";
827
828 function "*"
829 (Left : Real'Base;
830 Right : Complex_Matrix) return Complex_Matrix
831 renames Instantiations."*";
832
833 function "*"
834 (Left : Complex_Matrix;
835 Right : Real'Base) return Complex_Matrix
836 renames Instantiations."*";
837
838 ---------
839 -- "+" --
840 ---------
841
842 function "+" (Right : Complex_Vector) return Complex_Vector
843 renames Instantiations."+";
844
845 function "+"
846 (Left : Complex_Vector;
847 Right : Complex_Vector) return Complex_Vector
848 renames Instantiations."+";
849
850 function "+"
851 (Left : Real_Vector;
852 Right : Complex_Vector) return Complex_Vector
853 renames Instantiations."+";
854
855 function "+"
856 (Left : Complex_Vector;
857 Right : Real_Vector) return Complex_Vector
858 renames Instantiations."+";
859
860 function "+" (Right : Complex_Matrix) return Complex_Matrix
861 renames Instantiations."+";
862
863 function "+"
864 (Left : Complex_Matrix;
865 Right : Complex_Matrix) return Complex_Matrix
866 renames Instantiations."+";
867
868 function "+"
869 (Left : Real_Matrix;
870 Right : Complex_Matrix) return Complex_Matrix
871 renames Instantiations."+";
872
873 function "+"
874 (Left : Complex_Matrix;
875 Right : Real_Matrix) return Complex_Matrix
876 renames Instantiations."+";
877
878 ---------
879 -- "-" --
880 ---------
881
882 function "-"
883 (Right : Complex_Vector) return Complex_Vector
884 renames Instantiations."-";
885
886 function "-"
887 (Left : Complex_Vector;
888 Right : Complex_Vector) return Complex_Vector
889 renames Instantiations."-";
890
891 function "-"
892 (Left : Real_Vector;
893 Right : Complex_Vector) return Complex_Vector
894 renames Instantiations."-";
895
896 function "-"
897 (Left : Complex_Vector;
898 Right : Real_Vector) return Complex_Vector
899 renames Instantiations."-";
900
901 function "-" (Right : Complex_Matrix) return Complex_Matrix
902 renames Instantiations."-";
903
904 function "-"
905 (Left : Complex_Matrix;
906 Right : Complex_Matrix) return Complex_Matrix
907 renames Instantiations."-";
908
909 function "-"
910 (Left : Real_Matrix;
911 Right : Complex_Matrix) return Complex_Matrix
912 renames Instantiations."-";
913
914 function "-"
915 (Left : Complex_Matrix;
916 Right : Real_Matrix) return Complex_Matrix
917 renames Instantiations."-";
918
919 ---------
920 -- "/" --
921 ---------
922
923 function "/"
924 (Left : Complex_Vector;
925 Right : Complex) return Complex_Vector
926 renames Instantiations."/";
927
928 function "/"
929 (Left : Complex_Vector;
930 Right : Real'Base) return Complex_Vector
931 renames Instantiations."/";
932
933 function "/"
934 (Left : Complex_Matrix;
935 Right : Complex) return Complex_Matrix
936 renames Instantiations."/";
937
938 function "/"
939 (Left : Complex_Matrix;
940 Right : Real'Base) return Complex_Matrix
941 renames Instantiations."/";
942
943 -----------
944 -- "abs" --
945 -----------
946
947 function "abs" (Right : Complex_Vector) return Real'Base
948 renames Instantiations."abs";
949
950 --------------
951 -- Argument --
952 --------------
953
954 function Argument (X : Complex_Vector) return Real_Vector
955 renames Instantiations.Argument;
956
957 function Argument
958 (X : Complex_Vector;
959 Cycle : Real'Base) return Real_Vector
960 renames Instantiations.Argument;
961
962 function Argument (X : Complex_Matrix) return Real_Matrix
963 renames Instantiations.Argument;
964
965 function Argument
966 (X : Complex_Matrix;
967 Cycle : Real'Base) return Real_Matrix
968 renames Instantiations.Argument;
969
970 ----------------------------
971 -- Compose_From_Cartesian --
972 ----------------------------
973
974 function Compose_From_Cartesian (Re : Real_Vector) return Complex_Vector
975 renames Instantiations.Compose_From_Cartesian;
976
977 function Compose_From_Cartesian
978 (Re : Real_Vector;
979 Im : Real_Vector) return Complex_Vector
980 renames Instantiations.Compose_From_Cartesian;
981
982 function Compose_From_Cartesian (Re : Real_Matrix) return Complex_Matrix
983 renames Instantiations.Compose_From_Cartesian;
984
985 function Compose_From_Cartesian
986 (Re : Real_Matrix;
987 Im : Real_Matrix) return Complex_Matrix
988 renames Instantiations.Compose_From_Cartesian;
989
990 ------------------------
991 -- Compose_From_Polar --
992 ------------------------
993
994 function Compose_From_Polar
995 (Modulus : Real_Vector;
996 Argument : Real_Vector) return Complex_Vector
997 renames Instantiations.Compose_From_Polar;
998
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;
1004
1005 function Compose_From_Polar
1006 (Modulus : Real_Matrix;
1007 Argument : Real_Matrix) return Complex_Matrix
1008 renames Instantiations.Compose_From_Polar;
1009
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;
1015
1016 ---------------
1017 -- Conjugate --
1018 ---------------
1019
1020 function Conjugate (X : Complex_Vector) return Complex_Vector
1021 renames Instantiations.Conjugate;
1022
1023 function Conjugate (X : Complex_Matrix) return Complex_Matrix
1024 renames Instantiations.Conjugate;
1025
1026 -----------------
1027 -- Determinant --
1028 -----------------
1029
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;
1038
1039 -----------------
1040 -- Eigensystem --
1041 -----------------
1042
1043 procedure Eigensystem
1044 (A : Complex_Matrix;
1045 Values : out Real_Vector;
1046 Vectors : out Complex_Matrix)
1047 is
1048 N : constant Natural := Length (A);
1049
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)
1054 --
1055 -- is equivalent to the (2 * N, 2 * N) real eigenvalue problem:
1056 -- [ A, B ] [ u ] = Lambda * [ u ]
1057 -- [ -B, A ] [ v ] [ v ]
1058 --
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.
1061
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.
1065
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);
1069
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;
1084
1085 Eigensystem (M, Vals, Vecs);
1086
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);
1092
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;
1104
1105 -----------------
1106 -- Eigenvalues --
1107 -----------------
1108
1109 function Eigenvalues (A : Complex_Matrix) return Real_Vector is
1110 -- See Eigensystem for a description of the algorithm
1111
1112 N : constant Natural := Length (A);
1113 R : Real_Vector (A'Range (1));
1114
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;
1131
1132 Vals := Eigenvalues (M);
1133
1134 for J in 1 .. N loop
1135 R (A'First (1) + (J - 1)) := Vals (2 * J);
1136 end loop;
1137
1138 return R;
1139 end Eigenvalues;
1140
1141 --------
1142 -- Im --
1143 --------
1144
1145 function Im (X : Complex_Vector) return Real_Vector
1146 renames Instantiations.Im;
1147
1148 function Im (X : Complex_Matrix) return Real_Matrix
1149 renames Instantiations.Im;
1150
1151 -------------
1152 -- Inverse --
1153 -------------
1154
1155 function Inverse (A : Complex_Matrix) return Complex_Matrix is
1156 (Solve (A, Unit_Matrix (Length (A))));
1157
1158 -------------
1159 -- Modulus --
1160 -------------
1161
1162 function Modulus (X : Complex_Vector) return Real_Vector
1163 renames Instantiations.Modulus;
1164
1165 function Modulus (X : Complex_Matrix) return Real_Matrix
1166 renames Instantiations.Modulus;
1167
1168 --------
1169 -- Re --
1170 --------
1171
1172 function Re (X : Complex_Vector) return Real_Vector
1173 renames Instantiations.Re;
1174
1175 function Re (X : Complex_Matrix) return Real_Matrix
1176 renames Instantiations.Re;
1177
1178 ------------
1179 -- Set_Im --
1180 ------------
1181
1182 procedure Set_Im
1183 (X : in out Complex_Matrix;
1184 Im : Real_Matrix)
1185 renames Instantiations.Set_Im;
1186
1187 procedure Set_Im
1188 (X : in out Complex_Vector;
1189 Im : Real_Vector)
1190 renames Instantiations.Set_Im;
1191
1192 ------------
1193 -- Set_Re --
1194 ------------
1195
1196 procedure Set_Re
1197 (X : in out Complex_Matrix;
1198 Re : Real_Matrix)
1199 renames Instantiations.Set_Re;
1200
1201 procedure Set_Re
1202 (X : in out Complex_Vector;
1203 Re : Real_Vector)
1204 renames Instantiations.Set_Re;
1205
1206 -----------
1207 -- Solve --
1208 -----------
1209
1210 function Solve
1211 (A : Complex_Matrix;
1212 X : Complex_Vector) return Complex_Vector
1213 renames Instantiations.Solve;
1214
1215 function Solve
1216 (A : Complex_Matrix;
1217 X : Complex_Matrix) return Complex_Matrix
1218 renames Instantiations.Solve;
1219
1220 ---------------
1221 -- Transpose --
1222 ---------------
1223
1224 function Transpose
1225 (X : Complex_Matrix) return Complex_Matrix
1226 is
1227 R : Complex_Matrix (X'Range (2), X'Range (1));
1228 begin
1229 Transpose (X, R);
1230 return R;
1231 end Transpose;
1232
1233 -----------------
1234 -- Unit_Matrix --
1235 -----------------
1236
1237 function Unit_Matrix
1238 (Order : Positive;
1239 First_1 : Integer := 1;
1240 First_2 : Integer := 1) return Complex_Matrix
1241 renames Instantiations.Unit_Matrix;
1242
1243 -----------------
1244 -- Unit_Vector --
1245 -----------------
1246
1247 function Unit_Vector
1248 (Index : Integer;
1249 Order : Positive;
1250 First : Integer := 1) return Complex_Vector
1251 renames Instantiations.Unit_Vector;
1252
1253 end Ada.Numerics.Generic_Complex_Arrays;