1 ;; Copyright 2005, 2006 by Barton Willis
3 ;; This is free software; you can redistribute it and/or
4 ;; modify it under the terms of the GNU General Public License,
5 ;; http://www.gnu.org/copyleft/gpl.html.
7 ;; This software has NO WARRANTY, not even the implied warranty of
8 ;; MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
12 (macsyma-module conjugate
)
14 ($put
'$conjugate
1 '$version
)
15 ;; Let's remove built-in symbols from list for user-defined properties.
16 (setq $props
(remove '$conjugate $props
))
18 (defprop $conjugate tex-postfix tex
)
19 (defprop $conjugate
("^\\star") texsym
)
20 (defprop $conjugate
160. tex-lbp
)
21 (defprop $conjugate simp-conjugate operators
)
25 #-gcl
(:load-toplevel
:execute
)
26 (let (($context
'$global
) (context '$global
))
27 (meval '(($declare
) $conjugate $complex
))
28 ;; Let's remove built-in symbols from list for user-defined properties.
29 (setq $props
(remove '$conjugate $props
))))
31 ;; When a function commutes with the conjugate, give the function the
32 ;; commutes-with-conjugate property. The log function commutes with
33 ;; the conjugate on all of C except on the negative real axis. Thus
34 ;; log does not get the commutes-with-conjugate property. Instead,
35 ;; log gets the conjugate-function property.
37 ;; What important functions have I missed?
39 ;; (1) Arithmetic operators
41 (setf (get 'mplus
'commutes-with-conjugate
) t
)
42 (setf (get 'mtimes
'commutes-with-conjugate
) t
)
43 ;(setf (get 'mnctimes 'commutes-with-conjugate) t) ;; generally I think users will want this
44 (setf (get '%signum
'commutes-with-conjugate
) t
) ;; x=/=0, conjugate(signum(x)) = conjugate(x/abs(x)) = signum(conjugate(x))
45 ;; Trig-like functions and other such functions
47 (setf (get '%cosh
'commutes-with-conjugate
) t
)
48 (setf (get '%sinh
'commutes-with-conjugate
) t
)
49 (setf (get '%tanh
'commutes-with-conjugate
) t
)
50 (setf (get '%sech
'commutes-with-conjugate
) t
)
51 (setf (get '%csch
'commutes-with-conjugate
) t
)
52 (setf (get '%coth
'commutes-with-conjugate
) t
)
53 (setf (get '%cos
'commutes-with-conjugate
) t
)
54 (setf (get '%sin
'commutes-with-conjugate
) t
)
55 (setf (get '%tan
'commutes-with-conjugate
) t
)
56 (setf (get '%sec
'commutes-with-conjugate
) t
)
57 (setf (get '%csc
'commutes-with-conjugate
) t
)
58 (setf (get '%cot
'commutes-with-conjugate
) t
)
59 (setf (get '$atan2
'commutes-with-conjugate
) t
)
61 (setf (get '%jacobi_cn
'commutes-with-conjugate
) t
)
62 (setf (get '%jacobi_sn
'commutes-with-conjugate
) t
)
63 (setf (get '%jacobi_dn
'commutes-with-conjugate
) t
)
65 (setf (get '%gamma
'commutes-with-conjugate
) t
)
66 (setf (get '$pochhammer
'commutes-with-conjugate
) t
)
70 (setf (get '$matrix
'commutes-with-conjugate
) t
)
71 (setf (get 'mlist
'commutes-with-conjugate
) t
)
72 (setf (get '$set
'commutes-with-conjugate
) t
)
76 (setf (get 'mequal
'commutes-with-conjugate
) t
)
77 (setf (get 'mnotequal
'commutes-with-conjugate
) t
)
78 (setf (get '%transpose
'commutes-with-conjugate
) t
)
82 (setf (get '$max
'commutes-with-conjugate
) t
)
83 (setf (get '$min
'commutes-with-conjugate
) t
)
85 ;; When a function has the conjugate-function property,
86 ;; use a non-generic function to conjugate it. Not done:
87 ;; conjugate-functions for all the inverse trigonometric
90 ;; Trig like and hypergeometric like functions
92 (setf (get '%log
'conjugate-function
) 'conjugate-log
)
93 (setf (get 'mexpt
'conjugate-function
) 'conjugate-mexpt
)
94 (setf (get '%asin
'conjugate-function
) 'conjugate-asin
)
95 (setf (get '%acos
'conjugate-function
) 'conjugate-acos
)
96 (setf (get '%atan
'conjugate-function
) 'conjugate-atan
)
97 (setf (get '%atanh
'conjugate-function
) 'conjugate-atanh
)
99 ;;(setf (get '$asec 'conjugate-function) 'conjugate-asec)
100 ;;(setf (get '$acsc 'conjugate-function) 'conjugate-acsc)
101 (setf (get '%bessel_j
'conjugate-function
) 'conjugate-bessel-j
)
102 (setf (get '%bessel_y
'conjugate-function
) 'conjugate-bessel-y
)
103 (setf (get '%bessel_i
'conjugate-function
) 'conjugate-bessel-i
)
104 (setf (get '%bessel_k
'conjugate-function
) 'conjugate-bessel-k
)
106 (setf (get '%hankel_1
'conjugate-function
) 'conjugate-hankel-1
)
107 (setf (get '%hankel_2
'conjugate-function
) 'conjugate-hankel-2
)
111 (setf (get '%sum
'conjugate-function
) 'conjugate-sum
)
112 (setf (get '%product
'conjugate-function
) 'conjugate-product
)
114 ;; Return true iff Maxima can prove that z is not on the
115 ;; negative real axis.
117 (defun off-negative-real-axisp (z)
118 (setq z
(trisplit z
)) ; split into real and imaginary
119 (or (eq t
(mnqp (cdr z
) 0)) ; y # 0
120 (eq t
(mgqp (car z
) 0)))) ; x >= 0
122 (defun on-negative-real-axisp (z)
123 (setq z
(trisplit z
))
124 (and (eq t
(meqp (cdr z
) 0))
125 (eq t
(mgrp 0 (car z
)))))
127 (defun in-domain-of-asin (z)
128 (setq z
(trisplit z
))
129 (let ((x (car z
)) (y (cdr z
)))
130 (or (eq t
(mgrp y
0))
134 (eq t
(mgrp 1 x
))))))
136 ;; Return conjugate(log(x)). Actually, x is a lisp list (x).
138 (defun conjugate-log (x)
140 (cond ((off-negative-real-axisp x
)
141 (take '(%log
) (take '($conjugate
) x
)))
142 ((on-negative-real-axisp x
)
143 (add (take '(%log
) (neg x
)) (mul -
1 '$%i
'$%pi
)))
144 (t `(($conjugate simp
) ((%log simp
) ,x
)))))
146 ;; Return conjugate(x^p), where e = (x, p). Suppose x isn't on the negative real axis.
147 ;; Then conjugate(x^p) == conjugate(exp(p * log(x))) == exp(conjugate(p) * conjugate(log(x)))
148 ;; == exp(conjugate(p) * log(conjugate(x)) = conjugate(x)^conjugate(p). Thus, when
149 ;; x is off the negative real axis, commute the conjugate with ^. Also if p is an integer
150 ;; ^ commutes with the conjugate.
152 (defun conjugate-mexpt (e)
153 (let ((x (first e
)) (p (second e
)))
154 (if (or (off-negative-real-axisp x
) ($featurep p
'$integer
))
155 (power (take '($conjugate
) x
) (take '($conjugate
) p
))
156 `(($conjugate simp
) ,(power x p
)))))
158 (defun conjugate-sum (e)
159 (take '(%sum
) (take '($conjugate
) (first e
)) (second e
) (third e
) (fourth e
)))
161 (defun conjugate-product (e)
162 (take '(%product
) (take '($conjugate
) (first e
)) (second e
) (third e
) (fourth e
)))
164 (defun conjugate-asin (x)
166 (if (in-domain-of-asin x
) (take '(%asin
) (take '($conjugate
) x
))
167 `(($conjugate simp
) ((%asin
) ,x
))))
169 (defun conjugate-acos (x)
171 (if (in-domain-of-asin x
) (take '(%acos
) (take '($conjugate
) x
))
172 `(($conjugate simp
) ((%acos
) ,x
))))
174 (defun conjugate-atan (x)
177 (setq xx
(mul '$%i x
))
178 (if (in-domain-of-asin xx
)
179 (take '(%atan
) (take '($conjugate
) x
))
180 `(($conjugate simp
) ((%atan
) ,x
)))))
182 ;; atanh and asin are entire on the same set; see A&S Fig. 4.4 and 4.7.
184 (defun conjugate-atanh (x)
186 (if (in-domain-of-asin x
) (take '(%atanh
) (take '($conjugate
) x
))
187 `(($conjugate simp
) ((%atanh
) ,x
))))
189 ;; Integer order Bessel functions are entire; thus they commute with the
190 ;; conjugate (Schwartz refection principle). But non-integer order Bessel
191 ;; functions are not analytic along the negative real axis. Notice that A&S
192 ;; 9.1.40 isn't correct -- it says that the real order Bessel functions
193 ;; commute with the conjugate. Not true.
195 (defun conjugate-bessel-j (z)
196 (let ((n (first z
)) (x (second z
)))
197 (if (or ($featurep n
'$integer
) (off-negative-real-axisp x
))
198 (take '(%bessel_j
) (take '($conjugate
) n
) (take '($conjugate
) x
))
199 `(($conjugate simp
) ((%bessel_j simp
) ,@z
)))))
201 (defun conjugate-bessel-y (z)
202 (let ((n (first z
)) (x (second z
)))
203 (if (off-negative-real-axisp x
)
204 (take '(%bessel_y
) (take '($conjugate
) n
) (take '($conjugate
) x
))
205 `(($conjugate simp
) ((%bessel_y simp
) ,@z
)))))
207 (defun conjugate-bessel-i (z)
208 (let ((n (first z
)) (x (second z
)))
209 (if (or ($featurep n
'$integer
) (off-negative-real-axisp x
))
210 (take '(%bessel_i
) (take '($conjugate
) n
) (take '($conjugate
) x
))
211 `(($conjugate simp
) ((%bessel_i simp
) ,@z
)))))
213 (defun conjugate-bessel-k (z)
214 (let ((n (first z
)) (x (second z
)))
215 (if (off-negative-real-axisp x
)
216 (take '(%bessel_k
) (take '($conjugate
) n
) (take '($conjugate
) x
))
217 `(($conjugate simp
) ((%bessel_k simp
) ,@z
)))))
219 (defun conjugate-hankel-1 (z)
220 (let ((n (first z
)) (x (second z
)))
221 (if (off-negative-real-axisp x
)
222 (take '(%hankel_2
) (take '($conjugate
) n
) (take '($conjugate
) x
))
223 `(($conjugate simp
) ((%hankel_1 simp
) ,@z
)))))
225 (defun conjugate-hankel-2 (z)
226 (let ((n (first z
)) (x (second z
)))
227 (if (off-negative-real-axisp x
)
228 (take '(%hankel_1
) (take '($conjugate
) n
) (take '($conjugate
) x
))
229 `(($conjugate simp
) ((%hankel_2 simp
) ,@z
)))))
231 ;; When a function maps "everything" into the reals, put real-valued on the
232 ;; property list of the function name. This duplicates some knowledge that
233 ;; $rectform has. So it goes. The functions floor and ceiling also aren't
234 ;; defined off the real-axis. I suppose these functions could be given the
235 ;; real-valued property.
237 (setf (get '%imagpart
'real-valued
) t
)
238 (setf (get 'mabs
'real-valued
) t
)
239 (setf (get '%realpart
'real-valued
) t
)
240 (setf (get '%carg
'real-valued
) t
)
242 ;; manifestly-real-p isn't a great name, but it's OK. Since (manifestly-real-p '$inf) --> true
243 ;; it might be called manifestly-extended-real-p. A nonscalar isn't real.
245 ;; There might be some advantage to requiring that the subscripts to a $subvarp
246 ;; all be real. Why? Well li[n] maps reals to reals when n is real, but li[n] does
247 ;; not map the reals to reals when n is nonreal.
249 (defun manifestly-real-p (e)
252 (not (manifestly-pure-imaginary-p e
))
253 (not (manifestly-complex-p e
))
254 (not (manifestly-nonreal-p e
))
258 (and ($subvarp e
) (manifestly-real-p ($op e
)))))))
260 (defun manifestly-pure-imaginary-p (e)
266 (and (symbolp e
) (kindp e
'$imaginary
) (not ($nonscalarp e
)))
267 (and ($subvarp e
) (manifestly-pure-imaginary-p ($op e
)))))
268 ;; For now, let's use $csign on constant expressions only; once $csign improves,
269 ;; the ban on nonconstant expressions can be removed
270 (and ($constantp e
) (eq '$imaginary
($csign e
))))))
272 ;; Don't use (kindp e '$complex)!
274 (defun manifestly-complex-p (e)
276 (or (and (symbolp e
) (decl-complexp e
) (not ($nonscalarp e
)))
278 (and ($subvarp e
) (manifestly-complex-p ($op e
))
279 (not ($nonscalarp e
))))))
281 (defun manifestly-nonreal-p (e)
282 (and (symbolp e
) (or (member e
`($und $ind $zeroa $zerob t nil
)) ($nonscalarp e
))))
284 ;; For a subscripted function, conjugate always returns the conjugate noun-form.
285 ;; This could be repaired. For now, we don't have a scheme for conjugate(li[m](x)).
287 ;; We could make commutes_with_conjugate and maps_to_reals features. But I
288 ;; doubt it would get much use.
290 (defun simp-conjugate (e f z
)
292 (setq e
(simpcheck (cadr e
) z
)) ; simp and disrep if necessary
293 (cond ((complexp e
) (conjugate e
)) ; never happens, but might someday.
294 ((manifestly-real-p e
) e
)
295 ((manifestly-pure-imaginary-p e
) (mul -
1 e
))
296 ((manifestly-nonreal-p e
) `(($conjugate simp
) ,e
))
297 (($mapatom e
) `(($conjugate simp
) ,e
))
298 ((op-equalp e
'$conjugate
) (car (margs e
)))
300 ((and (symbolp (mop e
)) (get (mop e
) 'real-valued
)) e
)
302 ((and (symbolp (mop e
)) (get (mop e
) 'commutes-with-conjugate
))
303 (simplify (cons (list (mop e
)) (mapcar #'(lambda (s) (take '($conjugate
) s
)) (margs e
)))))
305 ((setq f
(and (symbolp (mop e
)) (get (mop e
) 'conjugate-function
)))
306 (funcall f
(margs e
)))
308 (t `(($conjugate simp
) ,e
))))
