Date: Thu, 18 May 89 11:24:50 EDT
From: "Mark Phillips" <mbp@lakisis.umd.edu>
Subject: Complex.m --- a package for calculations with complex expressions

Hello MathGroupers, 
     Guess what!?  This is not a bug report or a request to be removed
from this list!  

[[ applause here ]]

     It's a package for dealing with complex expressions.  Have you
ever noticed that Mathematica doesn't evaluate Re[x], Conjugate[x],
etc, unless x is a number?  There are many situations, however, in
which it is nice to be able to delcare a symbol to represent a real
quantity, so that complex expressions involving that symbol and Re, Im,
Conjugate, and so on can be simplified.  This package defines the
function DeclareReal which does just this.  It also defines several
other functions which are useful in the context of calculations with
complex expressions.

     I would be interested in hearing what people think of it, plus
suggestions for how to improve it.  One thing that I have not been able
to work out satisfactorily is how to tell Mathematica that Sqrt[x] is
real when x is real an nonnegative.  The package currently just checks
for N[x]>=0, which works if x is a real number but not otherwise.
Similar problems exist for other functions, most notably Power, of
which Sqrt is just a special case.  Any ideas?

Mark Phillips			  mbp@lakisis.umd.edu   (arpanet)
Department of Mathematics         (301) 454-2693
University of Maryland	          
College Park, Maryland  20742     

---------snip snip------snip snip------snip snip------snip snip---------
(*
 * Complex.m:
 *	Allows symbols to be declared to represent real numbers.
 *	Expressions involving these symbols and Re,Im,Conjugate, etc,
 *	then simplify nicely.  Also defines some useful functions for
 *	manipulating complex expressions.
 *
 * Author: Mark B. Phillips	(mbp@lakisis.umd.edu)
 *	   Dept. of Mathematics
 *	   University of Maryland
 *	   May 18, 1989
 *
 * Typical usage:      
 *      
 *      In[1]:= <<Complex.m
 *      
 *      In[2]:= DeclareReal[x,y]
 *      
 *      In[3]:= z = x + I y
 *      
 *      Out[3]= x + I y
 *      
 *      In[4]:= ComplexNormSquared[z]
 *      
 *               2    2
 *      Out[4]= x  + y
 *      
 *      In[5]:= z / Conjugate[z]
 *      
 *              x + I y
 *      Out[5]= --------
 *              x + -I y
 *      
 *      In[6]:= ComplexCanonify[%]
 *      
 *                         2    2
 *              2 I x y   x  - y
 *      Out[6]= ------- + -------
 *               2    2    2    2
 *              x  + y    x  + y
 *)

BeginPackage["Complex`"]

RealQ::usage =
"RealQ[x] yields true if x is a real expression or the name of a function
which takes real numbers to real numbers."

DeclareReal::usage =
"DeclareReal[x, y, ...] declares x, y, ... to be real variables or
functions which take real numbers to real numbers.  The x, y,...
must be atoms."

RealValuedQ::usage =
"RealValuedQ[f] yields true if f is the name of a real-valued function."

DeclareRealValuedQ::usage =
"DeclareRealValuedQ[f, g, ...] declares f, g, ... to be the names of
real-valued functions.  The f, g, ... must be atoms."

ComplexNorm::usage =
"ComplexNorm[z] gives the Euclidean norm of the complex quantity z.
It is mathematically equal to Abs[z], but is immediately evaluated
using the formula Sqrt[ComplexNormSquared[z]]."

ComplexNormSquared::usage =
"ComplexNormSquared[z] gives the square of the Euclidean norm of the
complex quantity z."

AbsReduce::usage =
"AbsReduce[expr] writes expr as Sqrt[ Re[expr]^2 + Im[expr]^2 ]."

ComplexParametricPlot::usage =
"ComplexParametricPlot[f, args] plots the real and imaginary
parts of f in a parametric plot; it is equivalent to
ParametricPlot[{Re[f],Im[f]}, args ]."

ComplexCanonify::usage =
"ComplexCanonify[expr] gives the complex expression expr in the form
a + I b, where a and b are real."

Begin["`private`"]

(*----------------------------------------------------------------------*)
(* This is ReIm.m as distributed with version 1.1 of Mathematica.  It   *)
(* is specifically included here (rather than included indirectly)      *)
(* because the 1.2 beta version works differently.                      *)

Unprotect[Re,Im]

RealQ[Pi] := True
RealQ[E] := True
RealQ[EulerGamma] := True
RealQ[Sin] := True
RealQ[Cos] := True
RealQ[Tan] := True
RealQ[n_Integer] := True
RealQ[n_Real] := True
RealQ[n_Rational] := True

Re[f_[x_]] := f[x] /; RealQ[f] && RealQ[x]

Re[x_] := x /; RealQ[x]
Im[x_] := 0 /; RealQ[x]

Re[x_+y_] := Re[x] + Re[y]
Im[x_+y_] := Im[x] + Im[y]

Re[x_ y_] := Re[x] Re[y] - Im[x] Im[y]
Im[x_ y_] := Re[x] Im[y] + Im[x] Re[y]

Re[ 1 / x_ ] := Re[x] / (Re[x]^2 + Im[x]^2)
Im[ 1 / x_ ] := -Im[x] / (Re[x]^2 + Im[x]^2)

Re[E^x_] := Cos[Im[x]] Exp[Re[x]]
Im[E^x_] := Sin[Im[x]] Exp[Re[x]]

Im[x_ ^ 2] := 2 Re[x] Im[x]

Re[ x_ ^ n_Integer ] := Block[{a, b},
	a = Round[n/2];
	b = n-a;
	Re[x^a] Re[x^b] - Im[x^a] Im[x^b]
	]

Im[ x_ ^ n_Integer ] := Block[{a, b},
	a = Round[n/2];
	b = n-a;
	Re[x^a] Im[x^b] + Im[x^a] Re[x^b]
	]

Re[x_Integer ^ n_Rational] := 0                /; IntegerQ[2n] && x<0
Im[x_Integer ^ n_Rational] := 
	(-x) ^ n   (-1) ^ ((Numerator[n]-1)/2) /; IntegerQ[2n] && x<0

Re[x_Integer ^ n_Rational] := x^n /; OddQ[Denominator[n]] || x>0
Im[x_Integer ^ n_Rational] := 0   /; OddQ[Denominator[n]] || x>0

Re[(-1) ^ n_Rational] := Cos[n Pi]
Im[(-1) ^ n_Rational] := Sin[n Pi]

Re[Sin[x_]] := Sin[Re[x]] Cosh[Im[x]]
Im[Sin[x_]] := Cos[Re[x]] Sinh[Im[x]]

Re[Cos[x_]] :=  Cos[Re[x]] Cosh[Im[x]]
Im[Cos[x_]] := -Sin[Re[x]] Sinh[Im[x]]

Re[Tan[x_]] := Re[Sin[x]/Cos[x]]
Im[Tan[x_]] := Im[Sin[x]/Cos[x]]

Protect[Re,Im]
(* End of ReIm.m                                                        *)
(*----------------------------------------------------------------------*)

RealQ[f_[x_]] := True /; RealValuedQ[f]
RealValuedQ[Re] := True
RealValuedQ[Im] := True
RealValuedQ[Abs] := True
RealValuedQ[Arg] := True

(* NOTE: RealQ should have an option to specify how it does
   comparisons for <, >, ==, etc. Right now we just use N[#]. *)
RealQ[Sqrt[x_]] := True /; (RealQ[x] && (N[x] >= 0))
RealQ[x_?RealQ + y_?RealQ] := True
RealQ[x_?RealQ * y_?RealQ] := True
RealQ[x_?RealQ ^ -1] := True

Unprotect[Conjugate]
Conjugate[x_] := x /; Re[x] == x
Conjugate[x_ + y_] := Conjugate[x] + Conjugate[y]
Conjugate[x_ * y_] := Conjugate[x] * Conjugate[y]
Conjugate[Conjugate[x_]] := x
Conjugate[Exp[y_]] := Exp[-y] /; Re[y] == 0
Protect[Conjugate]

Unprotect[Arg]
Arg[x_] := 0 /; (Im[x]==0 && N[x]>0)
Arg[x_] := Pi /; (Im[x]==0 && N[x]<0)
Arg[0] = 0
Protect[Arg]

SetAttributes[DeclareReal, Listable]
DeclareReal[x_,y__] := (DeclareReal[{x,y}] ; Null)
DeclareReal[x_] := (x/: RealQ[x] = True ; Null)

SetAttributes[DeclareRealValued, Listable]
DeclareRealValued[f_,g__] := (DeclareRealValued[{f,g}] ; Null)
DeclareRealValued[f_] := (f/: RealValuedQ[f] = True ; Null)

ComplexNormSquared[z_] := Expand[ z Conjugate[z] ]

ComplexNorm[z_] := Sqrt[ComplexNormSquared[z]]

AbsReduce[e_] := e //. Abs[z_] -> Sqrt[ Re[z]^2 + Im[z]^2 ]

ComplexParametricPlot[f_, args__] := ParametricPlot[ {Re[f],Im[f]}, args ]

ComplexCanonify[z_] := Together[Expand[Re[z]]] + I Together[Expand[Im[z]]]

End[]

EndPackage[]

