From: jzakiya <jzakiya@...>
Date: 2009-12-22T05:05:12+09:00
Subject: Re: Correcting complex math

On Dec 20, 1:05 pm, jzakiya <jzak...@gmail.com> wrote:
> On Dec 20, 9:51 am, Colin Bartlett <colin...@googlemail.com> wrote:
>
> > If this can be made to work, it might be useful.
> > Two thoughts that occur to me are:
>
> > 1. The code you show would probably need to use Rational(?),
> >    otherwise 1/n (or 1.0/n) may not be what one thinks it is.
>
> From a mathematical/engineering perspective it is
> understood that a root n is a positive integer > 1.
> In other words, what number 'r' multiplied by itself
> n times equals 'x'.
>
> A coding implementation can do the requisite checks on
> all the inputs and raise needed exceptions/warnings.
>
> Of course from an arithmetic perspective (a)^(1/n) can
> take any value of n (integer, float, complex) but if
> n not a positive integer (or interizable) number > 1
> then the discussion here may not apply to those cases.
>
> I think the syntax of (a)^(n^-1) vs (a)^(1/n) makes the
> point more mathematical clearer that what you are doing
> is finding the nth root of 'a' vs raising 'a' to a
> fractional power.
>
> > 2. What would one do with a number (Float or Rational)
> >    which is arbitrarily close to 1/n but which is not 1/n?
>
> If a number can be expressed as rational n = b/c then 'c' would be the
> root, and 'b' would be the power and
> (a)^b^(c^-1), otherwise see comment above.
>
>
>
> > I'd also like to take the opportunity of this thread
> > to suggest removing (!) something from Complex,
> > unless anyone can make a good case for retaining it:
>
> > require "complex"
> > n = -2                #=> -2
> > p = 1                 #=> 1
> > cn = Complex( n, 0 )  #=> Complex(-2, 0)
> > cp = Complex( p, 0 )  #=> Complex(1, 0)
> > n > p                 #=> false
> > cn > p                #=> true
> > cn > cp               #=> true
>
> > The Complex numbers are not an "ordered field",
> > at least not in any normally useful sense.
>
> > At the moment if either x or y is Complex
> > then  x <=> y  is calculated as  x.abs <=> y.abs;
>
> > # File complex.rb, line 312
> >   def <=> (other)
> >     self.abs <=> other.abs
> >   end
>
> > I think that can be misleading (see the examples above),
> > and I would prefer it if "<=>" and "<", etc,
> > were *not* defined for Complex,
> > and if one wants to compare absolute values
> > then one should do so explicitly.
>
> I would agree that <, >, are ambiguous. It's not just
> a matter of magnitude size but also angular 'size',
>
> Remember, complex variables are really vectors with
> magnitude and direction (angle).
>
> Take the example of a ferris wheel ride. Where you
> get on it is angular distance zero from your start.
> Once the ride starts, it will eventually be rotating at
> some constant speed, with changing angular direction.
>
> However, the max angular distance from the start with be
> when you are at the top of the ferris wheel, which is an
> angular distance of PI (180) degrees from your start. All
> other locations have CCP angular distance from the start.
>
> Thus, any complex number is defined with polar rep of
> |mag(a)| /_ angle(a)
>
> Therefore, to find all the n roots of real values, I can
> change my original code to make it more efficient and
> flexible by making it just find the magnitude and angles
> of the root locations.
>
> def Complex.root(a,n,k)
>   angle = (2*k+1)*PI/n    # for 'a' real negative
>   angle += PI/n if a > 0  # for 'a' real positive
>   mag = a.abs**n**-1
>   # The kth root of (a)**n**-1
>   mag*Complex(cos(angle),sin(angle))
> end
>
> The key is to find the appropriate angle rotation for
> 'a' positive or negative. Now by finding that angle
> and mag values, all the numeric classes can be written
> to find all the roots of real values with easy syntax.
>
> In Integer/Float you could now do:
>
> a.root(n) # default to +/- real root or first CCP
> a.root(n,k) # kth ccw root
> a.roots(n) # array of all (ccw) roots
>
> To make it intuitive (and follow Matz's POLS) let k
> roots go from 1st to kth, k = 1..n, which means that
> the routine that uses k just needs to subtract 1 from
> it before computing the angle.
>
> These features can be integrated within all the other
> numeric classes/modules as appropriate.

Example implementation for class Integer.
Output can be prettyfied if desired, but it's a start.

class Integer
  require 'complex'
  include Math

  def root(n,k=0) # k=1..n, k=0 default root
    return 'Invalid root: n not > 0' if n < 1
    return self if n == 1
    return 0 if self == 0
    mag = self.abs**n**-1
    if k==0  # return default roots
      return mag if self > 0
      return -1*mag if (n & 1) == 1 # n odd
      return mag*Complex(cos(PI/n),sin(PI/n))
    end

    # return kth root if root value
    k -= 1 # translate k = 1..n to k = 0..(n-1)
    angle = (2*k+1)*PI/n       # self negative
    angle += PI/n if self > 0  # self positive
    mag*Complex(cos(angle),sin(angle))
  end

  def roots(n)
    return 'Invalid root: n not > 0' if n < 1
    return self if n == 1
    return 0 if self == 0
    mag = self.abs**n**-1
    roots = [] ; theta = PI/n
    n.times do |k|
      angle=(2*k+1)*theta         # for neg integer
      angle += theta if self > 0  # for pos integer
      roots << mag*Complex(cos(angle), sin(angle))
    end
    return roots
  end
end