1module specfn;  % Special functions package for REDUCE.
2
3% Author:  Chris Cannam, Sept-Nov 1992.
4%          Winfried Neun, Nov 1992 ...
5%          contribution from various authors ...
6
7% Redistribution and use in source and binary forms, with or without
8% modification, are permitted provided that the following conditions are met:
9%
10%    * Redistributions of source code must retain the relevant copyright
11%      notice, this list of conditions and the following disclaimer.
12%    * Redistributions in binary form must reproduce the above copyright
13%      notice, this list of conditions and the following disclaimer in the
14%      documentation and/or other materials provided with the distribution.
15%
16% THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
17% AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO,
18% THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
19% PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNERS OR
20% CONTRIBUTORS
21% BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
22% CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
23% SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
24% INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
25% CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
26% ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
27% POSSIBILITY OF SUCH DAMAGE.
28%
29
30
31%  ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||  %
32%                                                                %
33%     Please report bugs to Winfried Neun,                       %
34%                           Konrad-Zuse-Zentrum                  %
35%                              fuer Informationstechnik Berlin,  %
36%                           Heilbronner Str. 10                  %
37%                           10711 Berlin - Wilmersdorf           %
38%                           Federal Republic of Germany          %
39%     or by email, neun@sc.ZIB-Berlin.de                         %
40%                                                                %
41%  ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||  %
42%                                                                %
43%     This package provides algebraic and numeric                %
44%     manipulations upon various special functions:              %
45%                                                                %
46%              -- Bernoulli Numbers                              %
47%              -- Gamma Function                                 %
48%              -- Pochhammer Notation                            %
49%              -- Digamma (Psi) Function and Derivatives         %
50%              -- Riemann Zeta Function                          %
51%              -- Bessel Functions J, Y, I and K                 %
52%              -- Airy Functions                                 %
53%              -- Hankel Functions H1 and H2                     %
54%              -- Kummer Hypergeometric Functions M and U        %
55%              -- Struve, Lommel and Whittaker Functions         %
56%              -- Integral funtions, Si, Ci, s_i (=si), Ei,...   %
57%              -- Simplification of Factorials                   %
58%              -- Solid and Spherical Harmonics                  %
59%              -- Jacobi Elliptic Functions                      %
60%              -- Elliptic Integrals                             %
61%              -- Jacobi Theta Functions                         %
62%              -- Weierstrassian Elliptic Functions              %
63%              -- Sigma Functions                                %
64%
65%     accessible through the new operators Bernoulli, Gamma,     %
66%     Pochhammer, Psi, Polygamma, Zeta, BesselJ, BesselY,        %
67%     BesselI, BesselK, Hankel1, Hankel2, KummerM, KummerU,      %
68%     AiryAi, AiryBi, AiryAiPrime, AiryBiPrime,                  %
69%     Jacobi{sn,cn,dn...}, Elliptic{E,F,K...}, JacobiE           %
70%     EllipticTheta{1,2,3,4}                                     %
71%     Weierstrass, WeierstrassZeta, sigma                        %
72%     sigma1, sigma2, sigma3                                     %
73%     Beta, StruveL, StruveH, Lommel1, Lommel2, WhittakerM       %
74%     and WhittakerW, with the new switch SaveSFs.               %
75%                                                                %
76%  ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||  %
77
78
79create!-package ('(specfn sfconsts sfgen sfbern dilog sfbinom sfpolys
80                   sfsums simpfact harmonic jsymbols recsimpl sfint
81		   sfellip sfellipi sftheta sfweier
82		  ),
83                 '(contrib specfn));
84
85exports sq2bf!*, c!:prec!:;
86
87switch savesfs=on;
88
89
90%symbolic inline procedure mksq!:new u;
91%  !*p2q(car fkern(u) .* 1);
92
93symbolic fluid '(bernoulli!-alist sf!-alist !*savesfs);
94
95symbolic ( bernoulli!-alist := nil );
96symbolic ( sf!-alist        := nil );
97
98symbolic inline procedure sq2bf!*(x);
99   (if fixp x then i2bf!: x
100      else ((if car y neq '!:rd!: then retag cdr !*rn2rd y
101               else retag cdr y) where y = !*a2f x));
102
103symbolic smacro procedure c!:prec!:;
104   !:bprec!:;
105
106
107% These functions are needed in other modules.
108%  complex!*on!*switch and complex!*off!*switch return t iff the
109%  switch complex was already in the correct position
110
111algebraic procedure complex!*on!*switch;
112  symbolic
113    if not !*complex then <<(onoff('complex,t) where !*msg := nil); nil>>
114     else t;
115
116algebraic procedure complex!*off!*switch;
117   if symbolic !*complex then
118      if symbolic !*msg then
119         << off msg; off complex; on msg >>
120      else off complex
121   else t;
122
123% complex!*restore!*switch takes the value returned by complex!*on!*switch or
124%  complex!*off!*switch and restore the switch complex to its former value,
125%  i.e. the switch is flipped if the argument is nil
126
127algebraic procedure complex!*restore!*switch(fl);
128  symbolic
129    begin scalar !*msg;
130      if not fl then onoff('complex,not !*complex)
131    end;
132
133%algebraic operator besselJ,besselY,besselI,besselK,hankel1,hankel2;
134%algebraic (operator kummerM, kummerU, struveh, struvel
135%                  ,lommel1, lommel2 ,whittakerm, whittakerw,
136%                   Airy_Ai, Airy_Bi,Airy_AiPrime,Airy_biprime);
137
138defautoload_operator(BesselJ,specbess);
139defautoload_operator(BesselY,specbess);
140defautoload_operator(BesselI,specbess);
141defautoload_operator(BesselK,specbess);
142defautoload_operator(hankel1,specbess);
143defautoload_operator(hankel2,specbess);
144defautoload_operator(KummerM,specbess);
145defautoload_operator(KummerU,specbess);
146defautoload_operator(StruveH,specbess);
147defautoload_operator(StruveL,specbess);
148defautoload_operator(lommel1,specbess);
149defautoload_operator(lommel2,specbess);
150defautoload_operator(WhittakerM,specbess);
151defautoload_operator(WhittakerW,specbess);
152defautoload_operator(Airy_Ai,specbess);
153defautoload_operator(Airy_Bi,specbess);
154defautoload_operator(Airy_Aiprime,specbess);
155defautoload_operator(Airy_Biprime,specbess);
156
157% elliptic functions and integrals
158%     defautoload_operator(jacobiam, specbess);
159%     defautoload_operator(jacobisn, specbess);
160%     defautoload_operator(jacobicn, specbess);
161%     defautoload_operator(jacobidn, specbess);
162%     defautoload_operator(jacobins, specbess);
163%     defautoload_operator(jacobinc, specbess);
164%     defautoload_operator(jacobind, specbess);
165%     defautoload_operator(jacobisc, specbess);
166%     defautoload_operator(jacobisd, specbess);
167%     defautoload_operator(jacobics, specbess);
168%     defautoload_operator(jacobids, specbess);
169%     defautoload_operator(jacobicd, specbess);
170%     defautoload_operator(jacobidc, specbess);
171%     defautoload_operator(jacobie,  specbess);
172%
173%     defautoload_operator(elliptice,       specbess);
174%     defautoload_operator(elliptice!',     specbess);
175%     defautoload_operator(ellipticf,       specbess);
176%     defautoload_operator(elliptick,       specbess);
177%     defautoload_operator(elliptick!',     specbess);
178%     defautoload_operator(elliptictheta1,  specbess);
179%     defautoload_operator(elliptictheta2,  specbess);
180%     defautoload_operator(elliptictheta3,  specbess);
181%     defautoload_operator(elliptictheta4,  specbess);
182
183%defautoload_operator(gamma,sfgamma);
184%defautoload_operator(igamma,sfgamma);
185%defautoload_operator(polygamma,sfgamma);
186%defautoload_operator(psi,sfgamma);
187%defautoload_operator(ibeta,sfgamma);
188%defautoload_operator(beta,sfgamma);
189%defautoload_operator(pochhammer,sfgamma);
190%defautoload_operator(zeta,sfgamma);
191
192endmodule;
193
194end;
195
196
197
198
199