p-adic field in ruby
Teichmüller character(sgn)とexpとlogも.
x != 0 に対して,
x = p^ord(x) sgn(x) e^log(x)
となるようにlogを拡張してる.
class Adic < Numeric def initialize adic, precision, arg = 0, ord = 0 @adic, @precision = adic, precision @modulus = adic ** precision @int, @ord = if arg.zero? [0, :infty] else while arg % adic == 0 arg /= adic ord += 1 end [arg % modulus, ord] end end attr_reader :adic, :precision, :modulus, :int, :ord def hash adic.hash ^ precision.hash ^ int.hash ^ ord.hash end def eql? b self.class.eql?(b.class) and adic.eql?(b.adic) and precision.eql?(b.precision) and int.eql?(b.int) and ord.eql?(b.ord) end def coerce b case b when self.class raise if adic != b.adic case when precision < b.precision [self.class.new(adic, precision, b.int, b.ord), self] when precision > b.precision [b, self.class.new(adic, b.precision, int, ord)] else [b, self] end when Integer [self.class.new(adic, precision, b), self] when Rational n, = coerce b.numerator d, = coerce b.denominator [n/d, self] else raise end end def == b case b when self.class, Rational, Integer b, a = coerce b a.eql? b else a, b = b.coerce self a == b end end def zero? int.zero? and ord == :infty end def zero self.class.new(adic, precision) end def unity self.class.new(adic, precision, 1) end def inspect x, s = int, [] precision.times do x, a = x.divmod adic s << a end b, e = if ord.eql?(:infty) ['', 0] elsif ord >= 0 ['0.1', ord] else ['10', -ord] end if adic <= 36 s[0].to_s(adic) + '.' + s[1..-1].map {|a| a.to_s(adic)}.join else s.inspect end + if e.zero? '' else ' * ' + b + if e == 1 '' else "^#{e}" end end end alias to_s inspect def + b case b when self.class, Rational, Integer b, a = coerce b case when a.zero? b when b.zero? a when b == - a a.zero else a, b = b, a unless b.ord >= a.ord self.class.new(adic, a.precision, a.int + b.int * adic ** (b.ord - a.ord), a.ord) end else a, b = b.coerce(self) a + b end end def -@ if zero? self else self.class.new(adic, precision, - int, ord) end end def - b self + (- b) end def * b case b when self.class, Integer, Rational b, a = coerce b if a.zero? or b.zero? a.zero else self.class.new(adic, a.precision, a.int * b.int, a.ord + b.ord) end else a, b = b.coerce(self) a * b end end def inverse raise ZeroDivisionError if zero? a, b, x, z = int, modulus, 1, 0 until b.zero? q, r = a.divmod b a, b, x, z = b, r, z, x - q * z end self.class.new(adic, precision, x, - ord) end def / b case b when self.class, Rational, Integer b, a = coerce b a * b.inverse else a, b = b.coerce(self) a / b end end def ** b if b < 0 b, s = -b, inverse else s = self end a = unity until b.zero? a *= s if b.odd? s *= s b >>= 1 end a end def unit raise Errno::EDOM.new 'unit' if zero? self.class.new(adic, precision, int) end def sgn if zero? zero else r = unit (precision - 1).times {r **= adic} r end end def exp 1 + exp0 end def log r = unit / sgn - (1 - r).log0 end def exp0 raise Errno::EDOM.new 'exp' unless zero? or (adic - 1) * ord > 1 r = self unless zero? n = 2 t = r * r / n while t.ord < precision + ord r += t n += 1 t *= self / n end end r end def log0 raise Errno::EDOM.new 'log' unless zero? or ord > 0 r = self unless zero? n = 2 t = r * r / n while t.ord < precision + ord r += t t *= n n += 1 t *= self / n end end r end end