p-adic-numbers/padic.py

152 lines
4.9 KiB
Python
Raw Normal View History

2017-05-14 16:07:53 +02:00
from fractions import Fraction
from sys import maxint
2017-05-14 16:27:28 +02:00
from modp import *
2017-05-14 16:07:53 +02:00
class PAdic:
2017-05-14 16:54:12 +02:00
def __init__(self, p):
self.p = p
2017-05-14 16:07:53 +02:00
self.val = '' # current known value
2017-05-15 09:48:33 +02:00
self.prec = 0 # current known precision, not containing trailing zeros
2017-05-14 16:07:53 +02:00
self.order = 0 # order/valuation of number
pass # initialize object as subclass perhaps
2017-05-15 09:48:33 +02:00
def get(self, prec, decimal = True):
2017-05-14 16:07:53 +02:00
'Return value of number with given precision.'
2017-05-14 16:27:28 +02:00
while self.prec < prec:
# update val based on value
2017-05-15 09:48:33 +02:00
self.val = str(int(self._nextdigit())) + self.val
2017-05-14 16:27:28 +02:00
self.prec += 1
2017-05-15 09:48:33 +02:00
if self.order < 0:
return (self.val[:self.order] + ('.' if decimal else '') + self.val[self.order:])[-prec-1:]
return (self.val + self.order * '0')[-prec:]
2017-05-14 16:27:28 +02:00
def _nextdigit(self):
'Calculate next digit of p-adic number.'
2017-05-14 16:07:53 +02:00
raise NotImplementedError
2017-05-14 16:54:12 +02:00
def getdigit(self, index):
'Return digit at given index.'
2017-05-15 09:48:33 +02:00
return int(self.get(index + 1, False)[0]) #int(self.get(index+1+int(index < -self.order))[-int(index < -self.order)])
2017-05-14 16:54:12 +02:00
2017-05-14 16:07:53 +02:00
# return value with precision up to 32 bits
def __int__(self):
2017-05-15 09:48:33 +02:00
return int(self.get(32), self.p)
2017-05-14 16:07:53 +02:00
def __str__(self):
return self.get(32)
# arithmetic operations
2017-05-15 09:48:33 +02:00
def __neg__(self):
return PAdicNeg(self.p, self)
2017-05-14 16:07:53 +02:00
def __add__(self, other):
2017-05-14 16:54:12 +02:00
return PAdicAdd(self.p, self, other)
2017-05-14 16:07:53 +02:00
def __radd__(self, other):
2017-05-14 16:54:12 +02:00
return PAdicAdd(self.p, other, self)
2017-05-14 16:07:53 +02:00
def __sub__(self, other):
2017-05-14 16:54:12 +02:00
return PAdicAdd(self.p, self, PAdicNeg(self.p, other))
2017-05-14 16:07:53 +02:00
def __rsub__(self, other):
2017-05-14 16:54:12 +02:00
return PAdicAdd(self.p, other, PAdicNeg(self.p, self))
2017-05-15 09:48:33 +02:00
def __mul__(self, other):
2017-05-14 16:54:12 +02:00
return PAdicMul(self.p, self, other)
2017-05-15 09:48:33 +02:00
def __rmul__(self, other):
2017-05-14 16:54:12 +02:00
return PAdicMul(self.p, other, self)
2017-05-14 16:07:53 +02:00
# p-adic norm
def __abs__(self):
if self.order == maxint:
return 0
2017-05-15 09:48:33 +02:00
numer = denom = 1
2017-05-14 16:07:53 +02:00
if self.order > 0:
2017-05-15 09:48:33 +02:00
numer = self.p ** self.order
2017-05-14 16:07:53 +02:00
if self.order < 0:
2017-05-15 09:48:33 +02:00
denom = self.p ** self.order
return Fraction(numer, denom)
2017-05-14 16:07:53 +02:00
class PAdicConst(PAdic):
2017-05-14 16:54:12 +02:00
def __init__(self, p, value):
PAdic.__init__(self, p)
2017-05-14 16:07:53 +02:00
value = Fraction(value)
# calculate valuation
if value == 0:
self.value = value
self.val = '0'
2017-05-14 18:09:52 +02:00
self.order = maxint
2017-05-14 16:07:53 +02:00
return
self.order = 0
2017-05-14 16:54:12 +02:00
while not value.numerator % self.p:
2017-05-14 16:07:53 +02:00
self.order += 1
2017-05-15 09:48:33 +02:00
value /= self.p
2017-05-14 16:54:12 +02:00
while not value.denominator % self.p:
2017-05-15 09:48:33 +02:00
self.order -= 1
value *= self.p
2017-05-14 16:07:53 +02:00
self.value = value
2017-05-15 09:48:33 +02:00
self.zero = not value
2017-05-14 16:07:53 +02:00
2017-05-15 09:48:33 +02:00
def get(self, prec, decimal = True):
if self.zero:
return '0' * prec
return PAdic.get(self, prec, decimal)
2017-05-14 18:09:52 +02:00
2017-05-14 16:27:28 +02:00
def _nextdigit(self):
'Calculate next digit of p-adic number.'
2017-05-14 16:54:12 +02:00
rem = ModP(self.p, self.value.numerator) / ModP(self.p, self.value.denominator)
2017-05-15 09:48:33 +02:00
self.value -= int(rem)
2017-05-14 16:54:12 +02:00
self.value /= self.p
return rem
class PAdicAdd(PAdic):
'Sum of two p-adic numbers.'
def __init__(self, p, arg1, arg2):
PAdic.__init__(self, p)
self.carry = 0
self.arg1 = arg1
self.arg2 = arg2
2017-05-15 09:48:33 +02:00
self.order = self.prec = min(arg1.order, arg2.order) # might be larger than this
arg1.order -= self.order
arg2.order -= self.order
2017-05-14 18:09:52 +02:00
# loop until first nonzero digit is found
2017-05-15 09:48:33 +02:00
self.index = 0
2017-05-14 18:09:52 +02:00
digit = self._nextdigit()
while digit == 0:
2017-05-15 09:48:33 +02:00
self.index += 1
2017-05-14 18:09:52 +02:00
self.order += 1
digit = self._nextdigit()
2017-05-15 09:48:33 +02:00
self.val += str(int(digit))
self.prec = 1
2017-05-14 16:54:12 +02:00
def _nextdigit(self):
2017-05-15 09:48:33 +02:00
s = self.arg1.getdigit(self.index) + self.arg2.getdigit(self.index) + self.carry
2017-05-14 16:54:12 +02:00
self.carry = s // self.p
2017-05-15 09:48:33 +02:00
self.index += 1
2017-05-14 16:54:12 +02:00
return s % self.p
2017-05-14 16:27:28 +02:00
2017-05-14 16:54:12 +02:00
class PAdicNeg(PAdic):
'Negation of a p-adic number.'
def __init__(self, p, arg):
PAdic.__init__(self, p)
self.arg = arg
2017-05-14 18:09:52 +02:00
self.order = arg.order
2017-05-14 16:54:12 +02:00
def _nextdigit(self):
if self.prec == 0:
2017-05-14 18:09:52 +02:00
return self.p - self.arg.getdigit(0) # cannot be p, 0th digit of arg must be nonzero
2017-05-15 09:48:33 +02:00
return self.p - 1 - self.arg.getdigit(self.prec)
2017-05-14 16:54:12 +02:00
class PAdicMul(PAdic):
'Product of two p-adic numbers.'
def __init__(self, p, arg1, arg2):
PAdic.__init__(self, p)
self.carry = 0
self.arg1 = arg1
self.arg2 = arg2
2017-05-15 09:48:33 +02:00
self.order = arg1.order + arg2.order
self.arg1.order = self.arg2.order = 0 # TODO requires copy
self.index = 0
2017-05-14 16:54:12 +02:00
def _nextdigit(self):
2017-05-15 09:48:33 +02:00
s = sum(self.arg1.getdigit(i) * self.arg2.getdigit(self.index - i) for i in xrange(self.index + 1)) + self.carry
2017-05-14 16:54:12 +02:00
self.carry = s // self.p
2017-05-15 09:48:33 +02:00
self.index += 1
2017-05-14 16:54:12 +02:00
return s % self.p