7 #include "cln/complex.h"
19 const cl_N log (const cl_N& a, const cl_N& b)
24 // (complex (/ (log (abs a)) (log b)) (/ (phase a) (log b))), genauer:
25 // falls a reell, >0: bekannt
26 // falls (= a 0): Error
27 // sonst: (phase a) errechnen, ein Float.
28 // b (falls rational) ins selbe Float-Format umwandeln,
29 // Imaginärteil := (/ (phase a) (log dieses_b)).
30 // Falls a rational: (log (abs a) b).
31 // Falls a komplex mit rationalem Real- und Imaginärteil,
32 // Betragsquadrat (expt (abs a) 2) exakt ausrechnen als
33 // (+ (expt (realpart a) 2) (expt (imagpart a) 2)).
34 // Setze Realteil := (/ (log Betragsquadrat b) 2).
35 // [Eventuell wird hierbei (log b) ein zweites Mal ausgerechnet,
36 // aber dies sowieso nur in Single-Precision.]
37 // Sonst bilde (abs a), ein Float, und (log (abs a)), ein Float,
38 // wandle b (falls rational) ins selbe Float-Format um,
39 // setze Realteil := (/ (log (abs a)) (log dieses_b)).
40 // sonst: (/ (log a) (log b))
48 // a und b sind beide reell und >0
51 // b ist reell und >0, a aber nicht.
53 // Imaginärteil (/ (phase a) (log b)) errechnen:
56 var cl_R angle = phase(a);
57 if (eq(angle,0)) // = Fixnum 0 <==> (= a 0) -> Error
58 { cl_error_division_by_0(); }
59 { DeclareType(cl_F,angle);
60 var cl_F bf = cl_somefloat(b,angle); // (float b)
64 // Realteil (/ (log (abs a)) (log b)) errechnen:
69 // a rational -> (log (abs a) b) errechnen:
70 re = log(abs(a),b); // NB: (abs a) > 0
75 if (rationalp(realpart(a)) && rationalp(imagpart(a))) {
76 // a komplex mit rationalem Real- und Imaginärteil a1,a2
77 var const cl_R& a1 = realpart(a);
78 var const cl_R& a2 = imagpart(a);
79 re = log(square(a1)+square(a2),b) / 2;
83 // Keine Chance für rationalen Realteil.
85 var cl_F abs_a = The(cl_F)(abs(a));
86 var cl_F log_abs_a = ln(abs_a);
87 var cl_F bf = cl_somefloat(b,log_abs_a); // (float b)
88 re = log_abs_a / ln(bf);
92 return complex_C(re,im);
96 // normaler komplexer Fall
97 return log(a) / log(b);