7 #include "cl_complex.h"
15 #include "cl_rational.h"
23 // Falls y=0: Ergebnis 1,
24 // [Nach CLTL folgendermaßen:
26 // x rational -> Fixnum 1
27 // x Float -> (float 1 x)
29 // x komplex rational -> Fixnum 1
30 // sonst: #C(1.0 0.0) im Float-Format des Real- bzw. Imaginärteils von x
32 // Falls x rational oder komplex rational oder |y| klein:
33 // x^|y| durch wiederholtes Quadrieren und Multiplizieren und evtl.
34 // Kehrwert-Bilden ermitteln.
35 // Sonst wie bei 'y Float'.
37 // Es gilt (expt x m/n) = (expt (expt x 1/n) m).
38 // Falls x in Q(i) liegt (also rational oder komplex rational ist):
39 // Sollte x^(m/n) in Q(i) liegen, so auch eine n-te Wurzel x^(1/n)
40 // (und bei n=2 oder n=4 damit auch alle n-ten Wurzeln x^(1/n) ).
41 // Falls x rational >=0: n-te Wurzel aus x nehmen. Ist sie rational,
42 // deren m-te Potenz als Ergebnis.
43 // Falls x rational <=0 oder komplex rational und n Zweierpotenz:
44 // n-te Wurzel aus x nehmen (mehrfaches sqrt). Ist sie rational oder
45 // komplex rational, deren m-te Potenz als Ergebnis.
46 // [Beliebige n betrachten!??]
47 // Falls n Zweierpotenz und |m|,n klein: n-te Wurzel aus x nehmen
48 // (mehrfaches sqrt), davon die m-te Potenz durch wiederholtes
49 // Quadrieren und Multiplizieren und evtl. Kehrwert-Bilden.
50 // Sonst wie bei 'y Float'.
51 // Falls y Float oder komplex:
53 // Falls Realteil von y >0 :
54 // liefere 0.0 falls x und y reell, #C(0.0 0.0) sonst.
57 // liefere 1.0 falls x und y reell, #C(1.0 0.0) sonst.
58 // Sonst: (exp (* (log x) y))
59 // Das Ergebnis liegt in Q(i), falls x in Q(i) liegt und 4y ein Integer ist.??
60 // Genauigkeit erhöhen, log2(|y|) Bits mehr??
61 // Bei x oder y rational und der andere Long-Float: bitte kein Single-Float!??
64 inline const cl_N expt_0 (const cl_N& x)
67 // y=0 -> 1 im Format von x.
79 var cl_R f = contagion(realpart(x),imagpart(x));
86 return complex_C(cl_float(1,f),cl_float(0,f));
90 // Exponent exakt 0 -> Ergebnis exakt 1
96 inline const cl_R contagion (const cl_N& x)
103 return contagion(realpart(x),imagpart(x));
107 const cl_N expt (const cl_N& x, const cl_N& y)
112 DeclareType(cl_RA,y);
118 return expt_0(x); // Liefere 1
119 if (fixnump(y)) // |y| klein ?
120 return expt(x,y); // exakt ausrechnen
124 DeclareType(cl_RA,x);
125 return expt(x,y); // exakt ausrechnen
129 if (rationalp(realpart(x)) && rationalp(imagpart(x)))
130 return expt(x,y); // exakt ausrechnen
132 // x nicht exakt und |y| groß
134 DeclareType(cl_RT,y);
136 var const cl_I& m = numerator(y);
137 var const cl_I& n = denominator(y);
141 DeclareType(cl_RA,x);
143 goto complex_rational;
146 if (rootp(x,n,&w)) // Wurzel rational?
151 if (rationalp(realpart(x)) && rationalp(imagpart(x)))
152 goto complex_rational;
157 var uintL k = power2p(n);
159 // n Zweierpotenz = 2^(k-1). n>1, also k>1
167 if (fixnump(m) && fixnump(n)) { // |m| und n klein?
168 var uintL _n = FN_to_UL(n);
169 if ((_n & (_n-1)) == 0) { // n Zweierpotenz?
171 until ((_n = _n >> 1) == 0)
179 // allgemeiner Fall (z.B. y Float oder komplex):
180 if (zerop(x)) { // x=0.0 ?
181 if (!plusp(realpart(y))) // Realteil von y <=0 ?
182 { cl_error_division_by_0(); }
183 if (realp(x) && realp(y)) {
186 var cl_R f = contagion(x,y);
187 // ein Float, da sonst x = Fixnum 0 gewesen wäre
188 { DeclareType(cl_F,f);
189 return cl_float(0,f); // 0.0
192 var cl_R f = contagion(contagion(x),contagion(y));
193 // ein Float, da sonst x = Fixnum 0 gewesen wäre
194 { DeclareType(cl_F,f);
195 var cl_F f0 = cl_float(0,f);
196 return complex_C(f0,f0); // #C(0.0 0.0)
200 if (zerop(y)) { // y=0.0 ?
201 if (realp(x) && realp(y)) {
204 var cl_R f = contagion(x,y);
205 // ein Float, da sonst y = Fixnum 0 gewesen wäre
206 { DeclareType(cl_F,f);
207 return cl_float(1,f);
210 var cl_R f = contagion(contagion(x),contagion(y));
211 // ein Float, da sonst y = Fixnum 0 gewesen wäre
212 { DeclareType(cl_F,f);
213 return complex_C(cl_float(1,f),cl_float(0,f));
217 return exp(log(x)*y);