3 # Complex numbers and associated mathematical functions
4 # -- Raphael Manfredi, Sept 1996
7 package Math::Complex; @ISA = qw(Exporter);
12 tan cotan asin acos atan acotan
13 sinh cosh tanh cotanh asinh acosh atanh acotanh
39 $package = 'Math::Complex'; # Package name
40 $display = 'cartesian'; # Default display format
43 # Object attributes (internal):
44 # cartesian [real, imaginary] -- cartesian form
45 # polar [rho, theta] -- polar form
46 # c_dirty cartesian form not up-to-date
47 # p_dirty polar form not up-to-date
48 # display display format (package's global when not set)
54 # Create a new complex number (cartesian form)
57 my $self = bless {}, shift;
59 $self->{cartesian} = [$re, $im];
68 # Create a new complex number (exponential form)
71 my $self = bless {}, shift;
72 my ($rho, $theta) = @_;
73 $theta += pi() if $rho < 0;
74 $self->{polar} = [abs($rho), $theta];
80 sub new { &make } # For backward compatibility only.
85 # Creates a complex number from a (re, im) tuple.
86 # This avoids the burden of writing Math::Complex->make(re, im).
90 return $package->make($re, $im);
96 # Creates a complex number from a (rho, theta) tuple.
97 # This avoids the burden of writing Math::Complex->emake(rho, theta).
100 my ($rho, $theta) = @_;
101 return $package->emake($rho, $theta);
107 # The number defined as 2 * pi = 360 degrees
110 $pi = 4 * atan2(1, 1) unless $pi;
117 # The number defined as i*i = -1;
120 $i = bless {} unless $i; # There can be only one i
121 $i->{cartesian} = [0, 1];
122 $i->{polar} = [1, pi/2];
129 # Attribute access/set routines
132 sub cartesian {$_[0]->{c_dirty} ? $_[0]->update_cartesian : $_[0]->{cartesian}}
133 sub polar {$_[0]->{p_dirty} ? $_[0]->update_polar : $_[0]->{polar}}
135 sub set_cartesian { $_[0]->{p_dirty}++; $_[0]->{cartesian} = $_[1] }
136 sub set_polar { $_[0]->{c_dirty}++; $_[0]->{polar} = $_[1] }
141 # Recompute and return the cartesian form, given accurate polar form.
143 sub update_cartesian {
145 my ($r, $t) = @{$self->{polar}};
146 $self->{c_dirty} = 0;
147 return $self->{cartesian} = [$r * cos $t, $r * sin $t];
154 # Recompute and return the polar form, given accurate cartesian form.
158 my ($x, $y) = @{$self->{cartesian}};
159 $self->{p_dirty} = 0;
160 return $self->{polar} = [0, 0] if $x == 0 && $y == 0;
161 return $self->{polar} = [sqrt($x*$x + $y*$y), atan2($y, $x)];
170 my ($z1, $z2, $regular) = @_;
171 my ($re1, $im1) = @{$z1->cartesian};
172 my ($re2, $im2) = ref $z2 ? @{$z2->cartesian} : ($z2);
173 unless (defined $regular) {
174 $z1->set_cartesian([$re1 + $re2, $im1 + $im2]);
177 return (ref $z1)->make($re1 + $re2, $im1 + $im2);
186 my ($z1, $z2, $inverted) = @_;
187 my ($re1, $im1) = @{$z1->cartesian};
188 my ($re2, $im2) = ref $z2 ? @{$z2->cartesian} : ($z2);
189 unless (defined $inverted) {
190 $z1->set_cartesian([$re1 - $re2, $im1 - $im2]);
194 (ref $z1)->make($re2 - $re1, $im2 - $im1) :
195 (ref $z1)->make($re1 - $re2, $im1 - $im2);
204 my ($z1, $z2, $regular) = @_;
205 my ($r1, $t1) = @{$z1->polar};
206 my ($r2, $t2) = ref $z2 ? @{$z2->polar} : (abs($z2), $z2 >= 0 ? 0 : pi);
207 unless (defined $regular) {
208 $z1->set_polar([$r1 * $r2, $t1 + $t2]);
211 return (ref $z1)->emake($r1 * $r2, $t1 + $t2);
220 my ($z1, $z2, $inverted) = @_;
221 my ($r1, $t1) = @{$z1->polar};
222 my ($r2, $t2) = ref $z2 ? @{$z2->polar} : (abs($z2), $z2 >= 0 ? 0 : pi);
223 unless (defined $inverted) {
224 $z1->set_polar([$r1 / $r2, $t1 - $t2]);
228 (ref $z1)->emake($r2 / $r1, $t2 - $t1) :
229 (ref $z1)->emake($r1 / $r2, $t1 - $t2);
235 # Computes z1**z2 = exp(z2 * log z1)).
238 my ($z1, $z2, $inverted) = @_;
239 return exp($z1 * log $z2) if defined $inverted && $inverted;
240 return exp($z2 * log $z1);
246 # Computes z1 <=> z2.
247 # Sorts on the real part first, then on the imaginary part. Thus 2-4i > 3+8i.
250 my ($z1, $z2, $inverted) = @_;
251 my ($re1, $im1) = @{$z1->cartesian};
252 my ($re2, $im2) = ref $z2 ? @{$z2->cartesian} : ($z2);
253 my $sgn = $inverted ? -1 : 1;
254 return $sgn * ($re1 <=> $re2) if $re1 != $re2;
255 return $sgn * ($im1 <=> $im2);
266 my ($r, $t) = @{$z->polar};
267 return (ref $z)->emake($r, pi + $t);
269 my ($re, $im) = @{$z->cartesian};
270 return (ref $z)->make(-$re, -$im);
276 # Compute complex's conjugate.
281 my ($r, $t) = @{$z->polar};
282 return (ref $z)->emake($r, -$t);
284 my ($re, $im) = @{$z->cartesian};
285 return (ref $z)->make($re, -$im);
291 # Compute complex's norm (rho).
295 my ($r, $t) = @{$z->polar};
302 # Compute complex's argument (theta).
306 return 0 unless ref $z;
307 my ($r, $t) = @{$z->polar};
314 # Compute sqrt(z) (positive only).
318 my ($r, $t) = @{$z->polar};
319 return (ref $z)->emake(sqrt($r), $t/2);
325 # Compute cbrt(z) (cubic root, primary only).
329 return $z ** (1/3) unless ref $z;
330 my ($r, $t) = @{$z->polar};
331 return (ref $z)->emake($r**(1/3), $t/3);
337 # Computes all nth root for z, returning an array whose size is n.
338 # `n' must be a positive integer.
340 # The roots are given by (for k = 0..n-1):
342 # z^(1/n) = r^(1/n) (cos ((t+2 k pi)/n) + i sin ((t+2 k pi)/n))
347 return undef unless $n > 0;
348 my ($r, $t) = ref $z ? @{$z->polar} : (abs($z), $z >= 0 ? 0 : pi);
351 my $theta_inc = 2 * pi / $n;
352 my $rho = $r ** (1/$n);
354 my $complex = ref($z) || $package;
355 for ($k = 0, $theta = $t / $n; $k < $n; $k++, $theta += $theta_inc) {
356 push(@root, $complex->emake($rho, $theta));
368 return $z unless ref $z;
369 my ($re, $im) = @{$z->cartesian};
380 return 0 unless ref $z;
381 my ($re, $im) = @{$z->cartesian};
392 my ($x, $y) = @{$z->cartesian};
393 return (ref $z)->emake(exp($x), $y);
403 my ($r, $t) = @{$z->polar};
404 return (ref $z)->make(log($r), $t);
414 $log10 = log(10) unless defined $log10;
415 return log($z) / $log10 unless ref $z;
416 my ($r, $t) = @{$z->polar};
417 return (ref $z)->make(log($r) / $log10, $t / $log10);
423 # Compute logn(z,n) = log(z) / log(n)
427 my $logn = $logn{$n};
428 $logn = $logn{$n} = log($n) unless defined $logn; # Cache log(n)
429 return log($z) / log($n);
435 # Compute cos(z) = (exp(iz) + exp(-iz))/2.
439 my ($x, $y) = @{$z->cartesian};
442 return (ref $z)->make(cos($x) * ($ey + $ey_1)/2, sin($x) * ($ey_1 - $ey)/2);
448 # Compute sin(z) = (exp(iz) - exp(-iz))/2.
452 my ($x, $y) = @{$z->cartesian};
455 return (ref $z)->make(sin($x) * ($ey + $ey_1)/2, cos($x) * ($ey - $ey_1)/2);
461 # Compute tan(z) = sin(z) / cos(z).
465 return sin($z) / cos($z);
471 # Computes cotan(z) = 1 / tan(z).
475 return cos($z) / sin($z);
481 # Computes the arc cosine acos(z) = -i log(z + sqrt(z*z-1)).
486 $cz = cplx($cz, 0) if !ref $cz && $cz < 0; # Force complex if <0
487 return ~i * log($z + sqrt $cz); # ~i is -i
493 # Computes the arc sine asin(z) = -i log(iz + sqrt(1-z*z)).
498 $cz = cplx($cz, 0) if !ref $cz && $cz < 0; # Force complex if <0
499 return ~i * log(i * $z + sqrt $cz); # ~i is -i
505 # Computes the arc tagent atan(z) = i/2 log((i+z) / (i-z)).
509 return i/2 * log((i + $z) / (i - $z));
515 # Computes the arc cotangent acotan(z) = -i/2 log((i+z) / (z-i))
519 return i/-2 * log((i + $z) / ($z - i));
525 # Computes the hyperbolic cosine cosh(z) = (exp(z) + exp(-z))/2.
529 my ($x, $y) = ref $z ? @{$z->cartesian} : ($z);
532 return ($ex + $ex_1)/2 unless ref $z;
533 return (ref $z)->make(cos($y) * ($ex + $ex_1)/2, sin($y) * ($ex - $ex_1)/2);
539 # Computes the hyperbolic sine sinh(z) = (exp(z) - exp(-z))/2.
543 my ($x, $y) = ref $z ? @{$z->cartesian} : ($z);
546 return ($ex - $ex_1)/2 unless ref $z;
547 return (ref $z)->make(cos($y) * ($ex - $ex_1)/2, sin($y) * ($ex + $ex_1)/2);
553 # Computes the hyperbolic tangent tanh(z) = sinh(z) / cosh(z).
557 return sinh($z) / cosh($z);
563 # Comptutes the hyperbolic cotangent cotanh(z) = cosh(z) / sinh(z).
567 return cosh($z) / sinh($z);
573 # Computes the arc hyperbolic cosine acosh(z) = log(z + sqrt(z*z-1)).
578 $cz = cplx($cz, 0) if !ref $cz && $cz < 0; # Force complex if <0
579 return log($z + sqrt $cz);
585 # Computes the arc hyperbolic sine asinh(z) = log(z + sqrt(z*z-1))
589 my $cz = $z*$z + 1; # Already complex if <0
590 return log($z + sqrt $cz);
596 # Computes the arc hyperbolic tangent atanh(z) = 1/2 log((1+z) / (1-z)).
600 my $cz = (1 + $z) / (1 - $z);
601 $cz = cplx($cz, 0) if !ref $cz && $cz < 0; # Force complex if <0
608 # Computes the arc hyperbolic cotangent acotanh(z) = 1/2 log((1+z) / (z-1)).
612 my $cz = (1 + $z) / ($z - 1);
613 $cz = cplx($cz, 0) if !ref $cz && $cz < 0; # Force complex if <0
620 # Compute atan(z1/z2).
623 my ($z1, $z2, $inverted) = @_;
624 my ($re1, $im1) = @{$z1->cartesian};
625 my ($re2, $im2) = ref $z2 ? @{$z2->cartesian} : ($z2);
627 if (defined $inverted && $inverted) { # atan(z2/z1)
628 return pi * ($re2 > 0 ? 1 : -1) if $re1 == 0 && $im1 == 0;
631 return pi * ($re1 > 0 ? 1 : -1) if $re2 == 0 && $im2 == 0;
641 # Set (fetch if no argument) display format for all complex numbers that
642 # don't happen to have overrriden it via ->display_format
644 # When called as a method, this actually sets the display format for
645 # the current object.
647 # Valid object formats are 'c' and 'p' for cartesian and polar. The first
648 # letter is used actually, so the type can be fully spelled out for clarity.
654 if (ref $self) { # Called as a method
656 } else { # Regular procedure call
662 return defined $self->{display} ? $self->{display} : $display
663 unless defined $format;
664 return $self->{display} = $format;
667 return $display unless defined $format;
668 return $display = $format;
674 # Show nicely formatted complex number under its cartesian or polar form,
675 # depending on the current display format:
677 # . If a specific display format has been recorded for this object, use it.
678 # . Otherwise, use the generic current default for all complex numbers,
679 # which is a package global variable.
686 $format = $z->{display} if defined $z->{display};
688 return $z->stringify_polar if $format =~ /^p/i;
689 return $z->stringify_cartesian;
693 # ->stringify_cartesian
695 # Stringify as a cartesian representation 'a+bi'.
697 sub stringify_cartesian {
699 my ($x, $y) = @{$z->cartesian};
702 $re = "$x" if abs($x) >= 1e-14;
703 if ($y == 1) { $im = 'i' }
704 elsif ($y == -1) { $im = '-i' }
705 elsif (abs($y) >= 1e-14) { $im = "${y}i" }
708 $str = $re if defined $re;
709 $str .= "+$im" if defined $im;
712 $str = '0' unless $str;
720 # Stringify as a polar representation '[r,t]'.
722 sub stringify_polar {
724 my ($r, $t) = @{$z->polar};
727 return '[0,0]' if $r <= 1e-14;
731 $nt = ($nt - int($nt)) * $tpi;
732 $nt += $tpi if $nt < 0; # Range [0, 2pi]
734 if (abs($nt) <= 1e-14) { $theta = 0 }
735 elsif (abs(pi-$nt) <= 1e-14) { $theta = 'pi' }
737 return "\[$r,$theta\]" if defined $theta;
740 # Okay, number is not a real. Try to identify pi/n and friends...
743 $nt -= $tpi if $nt > pi;
746 for ($k = 1, $kpi = pi; $k < 10; $k++, $kpi += pi) {
747 $n = int($kpi / $nt + ($nt > 0 ? 1 : -1) * 0.5);
748 if (abs($kpi/$n - $nt) <= 1e-14) {
749 $theta = ($nt < 0 ? '-':'').($k == 1 ? 'pi':"${k}pi").'/'.abs($n);
754 $theta = $nt unless defined $theta;
756 return "\[$r,$theta\]";
764 Math::Complex - complex numbers and associated mathematical functions
769 $z = Math::Complex->make(5, 6);
771 $j = cplxe(1, 2*pi/3);
775 This package lets you create and manipulate complex numbers. By default,
776 I<Perl> limits itself to real numbers, but an extra C<use> statement brings
777 full complex support, along with a full set of mathematical functions
778 typically associated with and/or extended to complex numbers.
780 If you wonder what complex numbers are, they were invented to be able to solve
781 the following equation:
785 and by definition, the solution is noted I<i> (engineers use I<j> instead since
786 I<i> usually denotes an intensity, but the name does not matter). The number
787 I<i> is a pure I<imaginary> number.
789 The arithmetics with pure imaginary numbers works just like you would expect
790 it with real numbers... you just have to remember that
796 5i + 7i = i * (5 + 7) = 12i
797 4i - 3i = i * (4 - 3) = i
802 Complex numbers are numbers that have both a real part and an imaginary
803 part, and are usually noted:
807 where C<a> is the I<real> part and C<b> is the I<imaginary> part. The
808 arithmetic with complex numbers is straightforward. You have to
809 keep track of the real and the imaginary parts, but otherwise the
810 rules used for real numbers just apply:
812 (4 + 3i) + (5 - 2i) = (4 + 5) + i(3 - 2) = 9 + i
813 (2 + i) * (4 - i) = 2*4 + 4i -2i -i*i = 8 + 2i + 1 = 9 + 2i
815 A graphical representation of complex numbers is possible in a plane
816 (also called the I<complex plane>, but it's really a 2D plane).
821 is the point whose coordinates are (a, b). Actually, it would
822 be the vector originating from (0, 0) to (a, b). It follows that the addition
823 of two complex numbers is a vectorial addition.
825 Since there is a bijection between a point in the 2D plane and a complex
826 number (i.e. the mapping is unique and reciprocal), a complex number
827 can also be uniquely identified with polar coordinates:
831 where C<rho> is the distance to the origin, and C<theta> the angle between
832 the vector and the I<x> axis. There is a notation for this using the
833 exponential form, which is:
837 where I<i> is the famous imaginary number introduced above. Conversion
838 between this form and the cartesian form C<a + bi> is immediate:
843 which is also expressed by this formula:
845 z = rho * exp(i * theta) = rho * (cos theta + i * sin theta)
847 In other words, it's the projection of the vector onto the I<x> and I<y>
848 axes. Mathematicians call I<rho> the I<norm> or I<modulus> and I<theta>
849 the I<argument> of the complex number. The I<norm> of C<z> will be
852 The polar notation (also known as the trigonometric
853 representation) is much more handy for performing multiplications and
854 divisions of complex numbers, whilst the cartesian notation is better
855 suited for additions and substractions. Real numbers are on the I<x>
856 axis, and therefore I<theta> is zero.
858 All the common operations that can be performed on a real number have
859 been defined to work on complex numbers as well, and are merely
860 I<extensions> of the operations defined on real numbers. This means
861 they keep their natural meaning when there is no imaginary part, provided
862 the number is within their definition set.
864 For instance, the C<sqrt> routine which computes the square root of
865 its argument is only defined for positive real numbers and yields a
866 positive real number (it is an application from B<R+> to B<R+>).
867 If we allow it to return a complex number, then it can be extended to
868 negative real numbers to become an application from B<R> to B<C> (the
869 set of complex numbers):
871 sqrt(x) = x >= 0 ? sqrt(x) : sqrt(-x)*i
873 It can also be extended to be an application from B<C> to B<C>,
874 whilst its restriction to B<R> behaves as defined above by using
875 the following definition:
877 sqrt(z = [r,t]) = sqrt(r) * exp(i * t/2)
879 Indeed, a negative real number can be noted C<[x,pi]>
880 (the modulus I<x> is always positive, so C<[x,pi]> is really C<-x>, a
882 and the above definition states that
884 sqrt([x,pi]) = sqrt(x) * exp(i*pi/2) = [sqrt(x),pi/2] = sqrt(x)*i
886 which is exactly what we had defined for negative real numbers above.
888 All the common mathematical functions defined on real numbers that
889 are extended to complex numbers share that same property of working
890 I<as usual> when the imaginary part is zero (otherwise, it would not
891 be called an extension, would it?).
893 A I<new> operation possible on a complex number that is
894 the identity for real numbers is called the I<conjugate>, and is noted
895 with an horizontal bar above the number, or C<~z> here.
902 z * ~z = (a + bi) * (a - bi) = a*a + b*b
904 We saw that the norm of C<z> was noted C<abs(z)> and was defined as the
905 distance to the origin, also known as:
907 rho = abs(z) = sqrt(a*a + b*b)
913 If z is a pure real number (i.e. C<b == 0>), then the above yields:
917 which is true (C<abs> has the regular meaning for real number, i.e. stands
918 for the absolute value). This example explains why the norm of C<z> is
919 noted C<abs(z)>: it extends the C<abs> function to complex numbers, yet
920 is the regular C<abs> we know when the complex number actually has no
921 imaginary part... This justifies I<a posteriori> our use of the C<abs>
922 notation for the norm.
926 Given the following notations:
928 z1 = a + bi = r1 * exp(i * t1)
929 z2 = c + di = r2 * exp(i * t2)
930 z = <any complex or real number>
932 the following (overloaded) operations are supported on complex numbers:
934 z1 + z2 = (a + c) + i(b + d)
935 z1 - z2 = (a - c) + i(b - d)
936 z1 * z2 = (r1 * r2) * exp(i * (t1 + t2))
937 z1 / z2 = (r1 / r2) * exp(i * (t1 - t2))
938 z1 ** z2 = exp(z2 * log z1)
940 abs(z1) = r1 = sqrt(a*a + b*b)
941 sqrt(z1) = sqrt(r1) * exp(i * t1/2)
942 exp(z1) = exp(a) * exp(i * b)
943 log(z1) = log(r1) + i*t1
944 sin(z1) = 1/2i (exp(i * z1) - exp(-i * z1))
945 cos(z1) = 1/2 (exp(i * z1) + exp(-i * z1))
947 atan2(z1, z2) = atan(z1/z2)
949 The following extra operations are supported on both real and complex
957 log10(z) = log(z) / log(10)
958 logn(z, n) = log(z) / log(n)
960 tan(z) = sin(z) / cos(z)
961 cotan(z) = 1 / tan(z)
963 asin(z) = -i * log(i*z + sqrt(1-z*z))
964 acos(z) = -i * log(z + sqrt(z*z-1))
965 atan(z) = i/2 * log((i+z) / (i-z))
966 acotan(z) = -i/2 * log((i+z) / (z-i))
968 sinh(z) = 1/2 (exp(z) - exp(-z))
969 cosh(z) = 1/2 (exp(z) + exp(-z))
970 tanh(z) = sinh(z) / cosh(z)
971 cotanh(z) = 1 / tanh(z)
973 asinh(z) = log(z + sqrt(z*z+1))
974 acosh(z) = log(z + sqrt(z*z-1))
975 atanh(z) = 1/2 * log((1+z) / (1-z))
976 acotanh(z) = 1/2 * log((1+z) / (z-1))
978 The I<root> function is available to compute all the I<n>th
979 roots of some complex, where I<n> is a strictly positive integer.
980 There are exactly I<n> such roots, returned as a list. Getting the
981 number mathematicians call C<j> such that:
985 is a simple matter of writing:
987 $j = ((root(1, 3))[1];
989 The I<k>th root for C<z = [r,t]> is given by:
991 (root(z, n))[k] = r**(1/n) * exp(i * (t + 2*k*pi)/n)
993 The I<spaceshift> operation is also defined. In order to ensure its
994 restriction to real numbers is conform to what you would expect, the
995 comparison is run on the real part of the complex number first,
996 and imaginary parts are compared only when the real parts match.
1000 To create a complex number, use either:
1002 $z = Math::Complex->make(3, 4);
1005 if you know the cartesian form of the number, or
1009 if you like. To create a number using the trigonometric form, use either:
1011 $z = Math::Complex->emake(5, pi/3);
1012 $x = cplxe(5, pi/3);
1014 instead. The first argument is the modulus, the second is the angle (in radians).
1015 (Mnmemonic: C<e> is used as a notation for complex numbers in the trigonometric
1018 It is possible to write:
1020 $x = cplxe(-3, pi/4);
1022 but that will be silently converted into C<[3,-3pi/4]>, since the modulus
1023 must be positive (it represents the distance to the origin in the complex
1026 =head1 STRINGIFICATION
1028 When printed, a complex number is usually shown under its cartesian
1029 form I<a+bi>, but there are legitimate cases where the polar format
1030 I<[r,t]> is more appropriate.
1032 By calling the routine C<Math::Complex::display_format> and supplying either
1033 C<"polar"> or C<"cartesian">, you override the default display format,
1034 which is C<"cartesian">. Not supplying any argument returns the current
1037 This default can be overridden on a per-number basis by calling the
1038 C<display_format> method instead. As before, not supplying any argument
1039 returns the current display format for this number. Otherwise whatever you
1040 specify will be the new display format for I<this> particular number.
1046 Math::Complex::display_format('polar');
1047 $j = ((root(1, 3))[1];
1048 print "j = $j\n"; # Prints "j = [1,2pi/3]
1049 $j->display_format('cartesian');
1050 print "j = $j\n"; # Prints "j = -0.5+0.866025403784439i"
1052 The polar format attempts to emphasize arguments like I<k*pi/n>
1053 (where I<n> is a positive integer and I<k> an integer within [-9,+9]).
1057 Thanks to overloading, the handling of arithmetics with complex numbers
1058 is simple and almost transparent.
1060 Here are some examples:
1064 $j = cplxe(1, 2*pi/3); # $j ** 3 == 1
1065 print "j = $j, j**3 = ", $j ** 3, "\n";
1066 print "1 + j + j**2 = ", 1 + $j + $j**2, "\n";
1068 $z = -16 + 0*i; # Force it to be a complex
1069 print "sqrt($z) = ", sqrt($z), "\n";
1071 $k = exp(i * 2*pi/3);
1072 print "$j - $k = ", $j - $k, "\n";
1076 Saying C<use Math::Complex;> exports many mathematical routines in the caller
1077 environment. This is construed as a feature by the Author, actually... ;-)
1079 The code is not optimized for speed, although we try to use the cartesian
1080 form for addition-like operators and the trigonometric form for all
1081 multiplication-like operators.
1083 The arg() routine does not ensure the angle is within the range [-pi,+pi]
1084 (a side effect caused by multiplication and division using the trigonometric
1087 All routines expect to be given real or complex numbers. Don't attempt to
1088 use BigFloat, since Perl has currently no rule to disambiguate a '+'
1089 operation (for instance) between two overloaded entities.
1093 Raphael Manfredi <F<Raphael_Manfredi@grenoble.hp.com>>