1------------------------------------------------------------------------------ 2-- -- 3-- GNAT COMPILER COMPONENTS -- 4-- -- 5-- ADA.NUMERICS.GENERIC_COMPLEX_ARRAYS -- 6-- -- 7-- B o d y -- 8-- -- 9-- Copyright (C) 2006-2018, 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------------------------------------------------------------------------------ 31 32with System.Generic_Array_Operations; use System.Generic_Array_Operations; 33 34package 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 First_1 => A'First (2), 1158 First_2 => A'First (1)))); 1159 1160 ------------- 1161 -- Modulus -- 1162 ------------- 1163 1164 function Modulus (X : Complex_Vector) return Real_Vector 1165 renames Instantiations.Modulus; 1166 1167 function Modulus (X : Complex_Matrix) return Real_Matrix 1168 renames Instantiations.Modulus; 1169 1170 -------- 1171 -- Re -- 1172 -------- 1173 1174 function Re (X : Complex_Vector) return Real_Vector 1175 renames Instantiations.Re; 1176 1177 function Re (X : Complex_Matrix) return Real_Matrix 1178 renames Instantiations.Re; 1179 1180 ------------ 1181 -- Set_Im -- 1182 ------------ 1183 1184 procedure Set_Im 1185 (X : in out Complex_Matrix; 1186 Im : Real_Matrix) 1187 renames Instantiations.Set_Im; 1188 1189 procedure Set_Im 1190 (X : in out Complex_Vector; 1191 Im : Real_Vector) 1192 renames Instantiations.Set_Im; 1193 1194 ------------ 1195 -- Set_Re -- 1196 ------------ 1197 1198 procedure Set_Re 1199 (X : in out Complex_Matrix; 1200 Re : Real_Matrix) 1201 renames Instantiations.Set_Re; 1202 1203 procedure Set_Re 1204 (X : in out Complex_Vector; 1205 Re : Real_Vector) 1206 renames Instantiations.Set_Re; 1207 1208 ----------- 1209 -- Solve -- 1210 ----------- 1211 1212 function Solve 1213 (A : Complex_Matrix; 1214 X : Complex_Vector) return Complex_Vector 1215 renames Instantiations.Solve; 1216 1217 function Solve 1218 (A : Complex_Matrix; 1219 X : Complex_Matrix) return Complex_Matrix 1220 renames Instantiations.Solve; 1221 1222 --------------- 1223 -- Transpose -- 1224 --------------- 1225 1226 function Transpose 1227 (X : Complex_Matrix) return Complex_Matrix 1228 is 1229 R : Complex_Matrix (X'Range (2), X'Range (1)); 1230 begin 1231 Transpose (X, R); 1232 return R; 1233 end Transpose; 1234 1235 ----------------- 1236 -- Unit_Matrix -- 1237 ----------------- 1238 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; 1244 1245 ----------------- 1246 -- Unit_Vector -- 1247 ----------------- 1248 1249 function Unit_Vector 1250 (Index : Integer; 1251 Order : Positive; 1252 First : Integer := 1) return Complex_Vector 1253 renames Instantiations.Unit_Vector; 1254 1255end Ada.Numerics.Generic_Complex_Arrays; 1256