123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159 |
- /* Magnitude of complex z
- z mag(z)
- - ------
- a a
- -a a
- (-1)^a 1
- exp(a + i b) exp(a)
- a b mag(a) mag(b)
- a + i b sqrt(a^2 + b^2)
- Notes
- 1. Handles mixed polar and rectangular forms, e.g. 1 + exp(i pi/3)
- 2. jean-francois.debroux reports that when z=(a+i*b)/(c+i*d) then
- mag(numerator(z)) / mag(denominator(z))
- must be used to get the correct answer. Now the operation is
- automatic.
- */
- #include "stdafx.h"
- #include "defs.h"
- void
- eval_mag(void)
- {
- push(cadr(p1));
- eval();
- mag();
- }
- void
- mag(void)
- {
- save();
- p1 = pop();
- push(p1);
- numerator();
- yymag();
- push(p1);
- denominator();
- yymag();
- divide();
- restore();
- }
- void
- yymag(void)
- {
- save();
- p1 = pop();
- if (isnegativenumber(p1)) {
- push(p1);
- negate();
- } else if (car(p1) == symbol(POWER) && equaln(cadr(p1), -1))
- // -1 to a power
- push_integer(1);
- else if (car(p1) == symbol(POWER) && cadr(p1) == symbol(E)) {
- // exponential
- push(caddr(p1));
- real();
- exponential();
- } else if (car(p1) == symbol(MULTIPLY)) {
- // product
- push_integer(1);
- p1 = cdr(p1);
- while (iscons(p1)) {
- push(car(p1));
- mag();
- multiply();
- p1 = cdr(p1);
- }
- } else if (car(p1) == symbol(ADD)) {
- // sum
- push(p1);
- rect(); // convert polar terms, if any
- p1 = pop();
- push(p1);
- real();
- push_integer(2);
- power();
- push(p1);
- imag();
- push_integer(2);
- power();
- add();
- push_rational(1, 2);
- power();
- simplify_trig();
- } else
- // default (all real)
- push(p1);
- restore();
- }
- #if SELFTEST
- static char *s[] = {
- "mag(a+i*b)",
- "(a^2+b^2)^(1/2)",
- "mag(exp(a+i*b))",
- "exp(a)",
- "mag(1)",
- "1",
- "mag(-1)",
- "1",
- "mag(1+exp(i*pi/3))",
- "3^(1/2)",
- "mag((a+i*b)/(c+i*d))",
- "(a^2+b^2)^(1/2)/((c^2+d^2)^(1/2))",
- "mag(exp(i theta))",
- "1",
- "mag(exp(-i theta))",
- "1",
- "mag((-1)^theta)",
- "1",
- "mag((-1)^(-theta))",
- "1",
- "mag(3*(-1)^theta)",
- "3",
- "mag(3*(-1)^(-theta))",
- "3",
- "mag(-3*(-1)^theta)",
- "3",
- "mag(-3*(-1)^(-theta))",
- "3",
- };
- void
- test_mag(void)
- {
- test(__FILE__, s, sizeof s / sizeof (char *));
- }
- #endif
|