arg.cpp 2.0 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138
  1. /* Argument (angle) of complex z
  2. z arg(z)
  3. - ------
  4. a 0
  5. -a -pi See note 3 below
  6. (-1)^a a pi
  7. exp(a + i b) b
  8. a b arg(a) + arg(b)
  9. a + i b arctan(b/a)
  10. Result by quadrant
  11. z arg(z)
  12. - ------
  13. 1 + i 1/4 pi
  14. 1 - i -1/4 pi
  15. -1 + i 3/4 pi
  16. -1 - i -3/4 pi
  17. Notes
  18. 1. Handles mixed polar and rectangular forms, e.g. 1 + exp(i pi/3)
  19. 2. Symbols in z are assumed to be positive and real.
  20. 3. Negative direction adds -pi to angle.
  21. Example: z = (-1)^(1/3), mag(z) = 1/3 pi, mag(-z) = -2/3 pi
  22. 4. jean-francois.debroux reports that when z=(a+i*b)/(c+i*d) then
  23. arg(numerator(z)) - arg(denominator(z))
  24. must be used to get the correct answer. Now the operation is
  25. automatic.
  26. */
  27. #include "stdafx.h"
  28. #include "defs.h"
  29. void
  30. eval_arg(void)
  31. {
  32. push(cadr(p1));
  33. eval();
  34. arg();
  35. }
  36. void
  37. arg(void)
  38. {
  39. save();
  40. p1 = pop();
  41. push(p1);
  42. numerator();
  43. yyarg();
  44. push(p1);
  45. denominator();
  46. yyarg();
  47. subtract();
  48. restore();
  49. }
  50. #define RE p2
  51. #define IM p3
  52. void
  53. yyarg(void)
  54. {
  55. save();
  56. p1 = pop();
  57. if (isnegativenumber(p1)) {
  58. push(symbol(PI));
  59. negate();
  60. } else if (car(p1) == symbol(POWER) && equaln(cadr(p1), -1)) {
  61. // -1 to a power
  62. push(symbol(PI));
  63. push(caddr(p1));
  64. multiply();
  65. } else if (car(p1) == symbol(POWER) && cadr(p1) == symbol(E)) {
  66. // exponential
  67. push(caddr(p1));
  68. imag();
  69. } else if (car(p1) == symbol(MULTIPLY)) {
  70. // product of factors
  71. push_integer(0);
  72. p1 = cdr(p1);
  73. while (iscons(p1)) {
  74. push(car(p1));
  75. arg();
  76. add();
  77. p1 = cdr(p1);
  78. }
  79. } else if (car(p1) == symbol(ADD)) {
  80. // sum of terms
  81. push(p1);
  82. rect();
  83. p1 = pop();
  84. push(p1);
  85. real();
  86. RE = pop();
  87. push(p1);
  88. imag();
  89. IM = pop();
  90. if (iszero(RE)) {
  91. push(symbol(PI));
  92. if (isnegative(IM))
  93. negate();
  94. } else {
  95. push(IM);
  96. push(RE);
  97. divide();
  98. arctan();
  99. if (isnegative(RE)) {
  100. push_symbol(PI);
  101. if (isnegative(IM))
  102. subtract(); // quadrant 1 -> 3
  103. else
  104. add(); // quadrant 4 -> 2
  105. }
  106. }
  107. } else
  108. // pure real
  109. push_integer(0);
  110. restore();
  111. }