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