P-Adic numbers, basic: Difference between revisions
m
Improbved code.
mNo edit summary |
m (Improbved code.) |
||
| (31 intermediate revisions by 8 users not shown) | |||
Line 50:
__TOC__
=={{header|C++}}==
This example displays p-adic numbers in standard mathematical format, consisting of a possibly infinite list of digits extending leftwards from the p-adic point. p-adic numbers are given corrrect to O(prime^40) and rational reconstructions are accurate to O(prime^20).
<syntaxhighlight lang="c++">
#include <cmath>
#include <cstdint>
#include <iostream>
#include <numeric>
#include <stdexcept>
#include <string>
#include <vector>
class Rational {
public:
Rational(const int32_t& aNumerator, const int32_t& aDenominator) {
if ( aDenominator < 0 ) {
numerator = -aNumerator;
denominator = -aDenominator;
} else {
numerator = aNumerator;
denominator = aDenominator;
}
if ( aNumerator == 0 ) {
denominator = 1;
}
const uint32_t divisor = std::gcd(numerator, denominator);
numerator /= divisor;
denominator /= divisor;
}
std::string to_string() const {
return std::to_string(numerator) + " / " + std::to_string(denominator);
}
private:
int32_t numerator;
int32_t denominator;
};
class P_adic {
public:
// Create a P_adic number, with p = 'prime', from the given rational 'numerator' / 'denominator'.
P_adic(const uint32_t& prime, int32_t numerator, int32_t denominator) : prime(prime) {
if ( denominator == 0 ) {
throw std::invalid_argument("Denominator cannot be zero");
}
order = 0;
// Process rational zero
if ( numerator == 0 ) {
digits.assign(DIGITS_SIZE, 0);
order = ORDER_MAX;
return;
}
// Remove multiples of 'prime' and adjust the order of the P_adic number accordingly
while ( modulo_prime(numerator) == 0 ) {
numerator /= static_cast<int32_t>(prime);
order += 1;
}
while ( modulo_prime(denominator) == 0 ) {
denominator /= static_cast<int32_t>(prime);
order -= 1;
}
// Standard calculation of P_adic digits
const uint64_t inverse = modulo_inverse(denominator);
while ( digits.size() < DIGITS_SIZE ) {
const uint32_t digit = modulo_prime(numerator * inverse);
digits.emplace_back(digit);
numerator -= digit * denominator;
if ( numerator != 0 ) {
// The denominator is not a power of a prime
uint32_t count = 0;
while ( modulo_prime(numerator) == 0 ) {
numerator /= static_cast<int32_t>(prime);
count += 1;
}
for ( uint32_t i = count; i > 1; --i ) {
digits.emplace_back(0);
}
}
}
}
// Return the sum of this P_adic number with the given P_adic number.
P_adic add(P_adic other) {
if ( prime != other.prime ) {
throw std::invalid_argument("Cannot add p-adic's with different primes");
}
std::vector<uint32_t> this_digits = digits;
std::vector<uint32_t> other_digits = other.digits;
std::vector<uint32_t> result;
// Adjust the digits so that the P_adic points are aligned
for ( int32_t i = 0; i < -order + other.order; ++i ) {
other_digits.insert(other_digits.begin(), 0);
}
for ( int32_t i = 0; i < -other.order + order; ++i ) {
this_digits.insert(this_digits.begin(), 0);
}
// Standard digit by digit addition
uint32_t carry = 0;
for ( uint32_t i = 0; i < std::min(this_digits.size(), other_digits.size()); ++i ) {
const uint32_t sum = this_digits[i] + other_digits[i] + carry;
const uint32_t remainder = sum % prime;
carry = ( sum >= prime ) ? 1 : 0;
result.emplace_back(remainder);
}
return P_adic(prime, result, all_zero_digits(result) ? ORDER_MAX : std::min(order, other.order));
}
// Return the Rational representation of this P_adic number.
Rational convert_to_rational() {
std::vector<uint32_t> numbers = digits;
// Zero
if ( numbers.empty() || all_zero_digits(numbers) ) {
return Rational(1, 0);
}
// Positive integer
if ( order >= 0 && ends_with(numbers, 0) ) {
for ( int32_t i = 0; i < order; ++i ) {
numbers.emplace(numbers.begin(), 0);
}
return Rational(convert_to_decimal(numbers), 1);
}
// Negative integer
if ( order >= 0 && ends_with(numbers, prime - 1) ) {
negate_digits(numbers);
for ( int32_t i = 0; i < order; ++i ) {
numbers.emplace(numbers.begin(), 0);
}
return Rational(-convert_to_decimal(numbers), 1);
}
// Rational
const P_adic copy(prime, digits, order);
P_adic sum(prime, digits, order);
int32_t denominator = 1;
do {
sum = sum.add(copy);
denominator += 1;
} while ( ! ( ends_with(sum.digits, 0) || ends_with(sum.digits, prime - 1) ) );
const bool negative = ends_with(sum.digits, 6);
if ( negative ) {
negate_digits(sum.digits);
}
int32_t numerator = negative ? -convert_to_decimal(sum.digits) : convert_to_decimal(sum.digits);
if ( order > 0 ) {
numerator *= std::pow(prime, order);
}
if ( order < 0 ) {
denominator *= std::pow(prime, -order);
}
return Rational(numerator, denominator);
}
// Return a string representation of this P_adic number.
std::string to_string() {
std::vector<uint32_t> numbers = digits;
pad_with_zeros(numbers);
std::string result = "";
for ( int64_t i = numbers.size() - 1; i >= 0; --i ) {
result += std::to_string(digits[i]);
}
if ( order >= 0 ) {
for ( int32_t i = 0; i < order; ++i ) {
result += "0";
}
result += ".0";
} else {
result.insert(result.length() + order, ".");
while ( result[result.length() - 1] == '0' ) {
result = result.substr(0, result.length() - 1);
}
}
return " ..." + result.substr(result.length() - PRECISION - 1);
}
private:
/**
* Create a P_adic, with p = 'prime', directly from a vector of digits.
*
* For example: with 'order' = 0, the vector [1, 2, 3, 4, 5] creates the p-adic ...54321.0,
* 'order' > 0 shifts the vector 'order' places to the left and
* 'order' < 0 shifts the vector 'order' places to the right.
*/
P_adic(const uint32_t& prime, const std::vector<uint32_t>& digits, const int32_t& order)
: prime(prime), digits(digits), order(order) {
}
// Transform the given vector of digits representing a P_adic number
// into a vector which represents the negation of the P_adic number.
void negate_digits(std::vector<uint32_t>& numbers) {
numbers[0] = modulo_prime(prime - numbers[0]);
for ( uint64_t i = 1; i < numbers.size(); ++i ) {
numbers[i] = prime - 1 - numbers[i];
}
}
// Return the multiplicative inverse of the given number modulo 'prime'.
uint32_t modulo_inverse(const uint32_t& number) const {
uint32_t inverse = 1;
while ( modulo_prime(inverse * number) != 1 ) {
inverse += 1;
}
return inverse;
}
// Return the given number modulo 'prime' in the range 0..'prime' - 1.
int32_t modulo_prime(const int64_t& number) const {
const int32_t div = static_cast<int32_t>(number % prime);
return ( div >= 0 ) ? div : div + prime;
}
// The given vector is padded on the right by zeros up to a maximum length of 'DIGITS_SIZE'.
void pad_with_zeros(std::vector<uint32_t>& vector) {
while ( vector.size() < DIGITS_SIZE ) {
vector.emplace_back(0);
}
}
// Return the given vector of base 'prime' integers converted to a decimal integer.
uint32_t convert_to_decimal(const std::vector<uint32_t>& numbers) const {
uint32_t decimal = 0;
uint32_t multiple = 1;
for ( const uint32_t& number : numbers ) {
decimal += number * multiple;
multiple *= prime;
}
return decimal;
}
// Return whether the given vector consists of all zeros.
bool all_zero_digits(const std::vector<uint32_t>& numbers) const {
for ( uint32_t number : numbers ) {
if ( number != 0 ) {
return false;
}
}
return true;
}
// Return whether the given vector ends with multiple instances of the given number.
bool ends_with(const std::vector<uint32_t>& numbers, const uint32_t& number) const {
for ( uint64_t i = numbers.size() - 1; i >= numbers.size() - PRECISION / 2; --i ) {
if ( numbers[i] != number ) {
return false;
}
}
return true;
}
uint32_t prime;
std::vector<uint32_t> digits;
int32_t order;
static const uint32_t PRECISION = 40;
static const uint32_t ORDER_MAX = 1'000;
static const uint32_t DIGITS_SIZE = PRECISION + 5;
};
int main() {
std::cout << "3-adic numbers:" << std::endl;
P_adic padic_one(3, -2, 87);
std::cout << "-2 / 87 => " << padic_one.to_string() << std::endl;
P_adic padic_two(3, 4, 97);
std::cout << "4 / 97 => " << padic_two.to_string() << std::endl;
P_adic sum = padic_one.add(padic_two);
std::cout << "sum => " << sum.to_string() << std::endl;
std::cout << "Rational = " << sum.convert_to_rational().to_string() << std::endl;
std::cout << std::endl;
std::cout << "7-adic numbers:" << std::endl;
padic_one = P_adic(7, 5, 8);
std::cout << "5 / 8 => " << padic_one.to_string() << std::endl;
padic_two = P_adic(7, 353, 30809);
std::cout << "353 / 30809 => " << padic_two.to_string() << std::endl;
sum = padic_one.add(padic_two);
std::cout << "sum => " << sum.to_string() << std::endl;
std::cout << "Rational = " << sum.convert_to_rational().to_string() << std::endl;
std::cout << std::endl;
}
</syntaxhighlight>
{{ out }}
<pre>
3-adic numbers:
-2 / 87 => ...101020111222001212021110002210102011122.2
4 / 97 => ...022220111100202001010001200002111122021.0
sum => ...201011000022210220101111202212220210220.2
Rational = 154 / 8439
7-adic numbers:
5 / 8 => ...424242424242424242424242424242424242425.0
353 / 30809 => ...560462505550343461155520004023663643455.0
sum => ...315035233123101033613062431266421216213.0
Rational = 156869 / 246472
</pre>
=={{header|FreeBASIC}}==
<
' ***********************************************
'subject: convert two rationals to p-adic numbers,
Line 413 ⟶ 739:
system
</syntaxhighlight>
{{out|Examples}}
<pre>
Line 600 ⟶ 926:
=={{header|Go}}==
{{trans|FreeBASIC}}
<
import "fmt"
Line 929 ⟶ 1,255:
fmt.Println()
}
}</
{{out}}
<pre>
2/1 + 0(2^4)
0 0 1 0
1/1 + 0(2^4)
0 0 0 1
+ =
0 0 1 1
3
4/1 + 0(2^4)
0 1 0 0
3/1 + 0(2^4)
0 0 1 1
+ =
0 1 1 1
-2/2
4/1 + 0(2^5)
0 0 1 0 0
3/1 + 0(2^5)
0 0 0 1 1
+ =
0 0 1 1 1
7
4/9 + 0(5^4)
4 2 1 1
8/9 + 0(5^4)
3 4 2 2
+ =
3 1 3 3
4/3
26/25 + 0(5^4)
0 1. 0 1
-109/125 + 0(5^4)
4. 0 3 1
+ =
0. 0 4 1
21/125
49/2 + 0(7^6)
3 3 3 4 0 0
-4851/2 + 0(7^6)
3 2 3 3 0 0
+ =
6 6 0 0 0 0
-2401
-9/5 + 0(3^8)
2 1 0 1 2 1 0 0
27/7 + 0(3^8)
1 2 0 1 1 0 0 0
+ =
1 0 1 0 0 1 0 0
72/35
5/19 + 0(2^12)
0 0 1 0 1 0 0 0 0 1 1 1
-101/384 + 0(2^12)
1 0 1 0 1. 0 0 0 1 0 0 1
+ =
1 1 1 0 0. 0 0 0 1 0 0 1
1/7296
2/7 + 0(10^7)
5 7 1 4 2 8 6
-1/7 + 0(10^7)
7 1 4 2 8 5 7
+ =
2 8 5 7 1 4 3
1/7
34/21 + 0(10^9)
9 5 2 3 8 0 9 5 4
-39034/791 + 0(10^9)
1 3 9 0 6 4 4 2 6
+ =
0 9 1 4 4 5 3 8 0
-16180/339
11/4 + 0(2^43)
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0. 1 1
679001/207 + 0(2^43)
0 1 0 0 0 1 0 1 0 1 0 0 0 0 0 1 1 0 0 0 1 0 1 1 1 1 0 0 0 0 0 1 0 1 0 0 1 0 1 0 1 1 1
+ =
0 0 0 1 0 1 0 1 0 0 0 0 0 1 1 0 0 0 1 0 1 1 1 1 0 0 0 0 0 1 0 1 0 0 1 0 1 1 0 0 1. 1 1
2718281/828
-8/9 + 0(23^9)
2 12 17 20 10 5 2 12 17
302113/92 + 0(23^9)
5 17 5 17 6 0 10 12. 2
+ =
18 12 3 4 11 3 0 6. 2
2718281/828
-22/7 + 0(3^23)
1 0 2 1 2 0 1 0 2 1 2 0 1 0 2 1 2 0 1 0 2 0 2
46071/379 + 0(3^23)
2 0 1 2 1 2 1 2 2 1 2 1 0 0 2 2 0 1 1 2 1 0 0
+ =
0 1 1 1 1 0 0 0 2 0 1 1 1 1 2 0 2 2 0 0 0 0 2
314159/2653
-22/7 + 0(32749^3)
28070 18713 23389
46071/379 + 0(32749^3)
4493 8727 10145
+ =
32563 27441 785
314159/2653
35/61 + 0(5^20)
2 3 2 3 0 2 4 1 3 3 0 0 4 0 2 2 1 2 2 0
9400/109 + 0(5^20)
3 1 4 4 1 2 3 4 4 3 4 1 1 3 1 1 2 4 0 0
+ =
1 0 2 2 2 0 3 1 3 1 4 2 0 3 3 3 4 1 2 0
577215/6649
-101/109 + 0(61^7)
33 1 7 16 48 7 50
583376/6649 + 0(61^7)
33 1 7 16 49 34. 35
+ =
34 8 24 3 57 23. 35
577215/6649
-25/26 + 0(7^13)
2 6 5 0 5 4 4 0 1 6 1 2 2
5571/137 + 0(7^13)
3 2 4 1 4 5 4 2 2 5 5 3 5
+ =
6 2 2 2 3 3 1 2 4 4 6 6 0
141421/3562
1/4 + 0(7^11)
1 5 1 5 1 5 1 5 1 5 2
9263/2837 + 0(7^11)
6 5 6 6 0 3 2 0 4 4 1
+ =
1 4 1 4 2 1 3 5 6 2 3
39889/11348
122/407 + 0(7^11)
6 2 0 3 0 6 2 4 4 4 3
-517/1477 + 0(7^11)
1 2 3 4 3 5 4 6 4 1. 1
+ =
3 2 6 5 3 1 2 4 1 4. 1
-27584/90671
5/8 + 0(7^11)
4 2 4 2 4 2 4 2 4 2 5
353/30809 + 0(7^11)
2 3 6 6 3 6 4 3 4 5 5
+ =
6 6 4 2 1 2 1 6 2 1 3
47099/10977
</pre>
=={{header|Haskell}}==
p-Adic numbers in this implementation are represented in floating point manner, with p-adic unit as mantissa and p-adic norm as an exponent. The base is encoded implicitly at type level, so that combination of p-adics with different bases won't typecheck.
Textual presentation is given in traditional form with infinite part going to the left.
p-Adic arithmetics and conversion between rationals is implemented as instances of <code>Eq</code>, <code>Num</code>, <code>Fractional</code> and <code>Real</code> classes, so, they could be treated as usual real numbers (up to existence of some rationals for non-prime bases).
<syntaxhighlight lang="haskell">{-# LANGUAGE KindSignatures, DataKinds #-}
module Padic where
import Data.Ratio
import Data.List (genericLength)
import GHC.TypeLits
data Padic (n :: Nat) = Null
| Padic { unit :: [Int], order :: Int }
-- valuation of the base
modulo :: (KnownNat p, Integral i) => Padic p -> i
modulo = fromIntegral . natVal
-- Constructor for zero value
pZero :: KnownNat p => Padic p
pZero = Padic (repeat 0) 0
-- Smart constructor, adjusts trailing zeros with the order.
mkPadic :: (KnownNat p, Integral i) => [i] -> Int -> Padic p
mkPadic u k = go 0 (fromIntegral <$> u)
where
go 17 _ = pZero
go i (0:u) = go (i+1) u
go i u = Padic u (k-i)
-- Constructor for p-adic unit
mkUnit :: (KnownNat p, Integral i) => [i] -> Padic p
mkUnit u = mkPadic u 0
-- Zero test (up to 1/p^17)
isZero :: KnownNat p => Padic p -> Bool
isZero (Padic u _) = all (== 0) (take 17 u)
isZero _ = False
-- p-adic norm
pNorm :: KnownNat p => Padic p -> Ratio Int
pNorm Null = undefined
pNorm p = fromIntegral (modulo p) ^^ (- order p)
-- test for an integerness up to p^-17
isInteger :: KnownNat p => Padic p -> Bool
isInteger Null = False
isInteger (Padic s k) = case splitAt k s of
([],i) -> length (takeWhile (==0) $ reverse (take 20 i)) > 3
_ -> False
-- p-adics are shown with 1/p^17 precision
instance KnownNat p => Show (Padic p) where
show Null = "Null"
show x@(Padic u k) =
show (modulo x) ++ "-adic: " ++
(case si of {[] -> "0"; _ -> si})
++ "." ++
(case f of {[] -> "0"; _ -> sf})
where
(f,i) = case compare k 0 of
LT -> ([], replicate (-k) 0 ++ u)
EQ -> ([], u)
GT -> splitAt k (u ++ repeat 0)
sf = foldMap showD $ reverse $ take 17 f
si = foldMap showD $ dropWhile (== 0) $ reverse $ take 17 i
el s = if length s > 16 then "…" else ""
showD n = [(['0'..'9']++['a'..'z']) !! n]
instance KnownNat p => Eq (Padic p) where
a == b = isZero (a - b)
instance KnownNat p => Ord (Padic p) where
compare = error "Ordering is undefined fo p-adics."
instance KnownNat p => Num (Padic p) where
fromInteger 0 = pZero
fromInteger n = pAdic (fromInteger n)
x@(Padic a ka) + Padic b kb = mkPadic s k
where
k = ka `max` kb
s = addMod (modulo x)
(replicate (k-ka) 0 ++ a)
(replicate (k-kb) 0 ++ b)
_ + _ = Null
x@(Padic a ka) * Padic b kb =
mkPadic (mulMod (modulo x) a b) (ka + kb)
_ * _ = Null
negate x@(Padic u k) =
case map (\y -> modulo x - 1 - y) u of
n:ns -> Padic ((n+1):ns) k
[] -> pZero
negate _ = Null
abs p = pAdic (pNorm p)
signum = undefined
------------------------------------------------------------
-- conversion from rationals to p-adics
instance KnownNat p => Fractional (Padic p) where
fromRational = pAdic
recip Null = Null
recip x@(Padic (u:us) k)
| isZero x = Null
| gcd p u /= 1 = Null
| otherwise = mkPadic res (-k)
where
p = modulo x
res = longDivMod p (1:repeat 0) (u:us)
pAdic :: (Show i, Integral i, KnownNat p)
=> Ratio i -> Padic p
pAdic 0 = pZero
pAdic x = res
where
p = modulo res
(k, q) = getUnit p x
(n, d) = (numerator q, denominator q)
res = maybe Null process $ recipMod p d
process r = mkPadic (series n) k
where
series n
| n == 0 = repeat 0
| n `mod` p == 0 = 0 : series (n `div` p)
| otherwise =
let m = (n * r) `mod` p
in m : series ((n - m * d) `div` p)
------------------------------------------------------------
-- conversion from p-adics to rationals
-- works for relatively small denominators
instance KnownNat p => Real (Padic p) where
toRational Null = error "no rational representation!"
toRational x@(Padic s k) = res
where
p = modulo x
res = case break isInteger $ take 10000 $ iterate (x +) x of
(_,[]) -> - toRational (- x)
(d, i:_) -> (fromBase p (unit i) * (p^(- order i))) % (genericLength d + 1)
fromBase p = foldr (\x r -> r*p + x) 0 .
take 20 . map fromIntegral
--------------------------------------------------------------------------------
-- helper functions
-- extracts p-adic unit from a rational number
getUnit :: Integral i => i -> Ratio i -> (Int, Ratio i)
getUnit p x = (genericLength k1 - genericLength k2, c)
where
(k1,b:_) = span (\n -> denominator n `mod` p == 0) $
iterate (* fromIntegral p) x
(k2,c:_) = span (\n -> numerator n `mod` p == 0) $
iterate (/ fromIntegral p) b
-- Reciprocal of a number modulo p (extended Euclidean algorithm).
-- For non-prime p returns Nothing non-invertible element of the ring.
recipMod :: Integral i => i -> i -> Maybe i
recipMod p 1 = Just 1
recipMod p a | gcd p a == 1 = Just $ go 0 1 p a
| otherwise = Nothing
where
go t _ _ 0 = t `mod` p
go t nt r nr =
let q = r `div` nr
in go nt (t - q*nt) nr (r - q*nr)
-- Addition of two sequences modulo p
addMod p = go 0
where
go 0 [] ys = ys
go 0 xs [] = xs
go s [] ys = go 0 [s] ys
go s xs [] = go 0 xs [s]
go s (x:xs) (y:ys) =
let (q, r) = (x + y + s) `divMod` p
in r : go q xs ys
-- Subtraction of two sequences modulo p
subMod p a (b:bs) = addMod p a $ (p-b) : ((p - 1 -) <$> bs)
-- Multiplication of two sequences modulo p
mulMod p as [b] = mulMod p [b] as
mulMod p as bs = case as of
[0] -> repeat 0
[1] -> bs
[a] -> go 0 bs
where
go s [] = [s]
go s (b:bs) =
let (q, r) = (a * b + s) `divMod` p
in r : go q bs
as -> go bs
where
go [] = []
go (b:bs) =
let c:cs = mulMod p [b] as
in c : addMod p (go bs) cs
-- Division of two sequences modulo p
longDivMod p a (b:bs) = case recipMod p b of
Nothing -> error $
show b ++ " is not invertible modulo " ++ show p
Just r -> go a
where
go [] = []
go (0:xs) = 0 : go xs
go (x:xs) =
let m = (x*r) `mod` p
_:zs = subMod p (x:xs) (mulMod p [m] (b:bs))
in m : go zs</syntaxhighlight>
Convertation between rationals and p-adic numbers
<pre>:set -XDataKinds
λ> 0 :: Padic 5
5-adic: 0.0
λ> 3 :: Padic 5
5-adic: 3.0
λ> 1/3 :: Padic 5
5-adic: 13131313131313132.0
λ> 0.08 :: Padic 5
5-adic: 0.02
λ> 0.08 :: Padic 2
2-adic: 1011100001010010.0
λ> 0.08 :: Padic 7
7-adic: 36303630363036304.0
λ> toRational it
2 % 25
λ> λ> 25 :: Padic 10
10-adic: 25.0
λ> 1/25 :: Padic 10
Null
λ> -12/23 :: Padic 3
3-adic: 21001011200210010.0
λ> toRational it
(-12) % 23</pre>
Arithmetic:
<pre>λ> pAdic (12/25) + pAdic (23/56) :: Padic 7
7-adic: 32523252325232524.2
λ> pAdic (12/25 + 23/56) :: Padic 7
7-adic: 32523252325232524.2
λ> toRational it
1247 % 1400
λ> 12/25 + 23/56 :: Rational
1247 % 1400
λ> let x = 2/7 :: Padic 13
λ> (2*x - x^2) / 3
13-adic: 35ab53c972179036.0λ
> toRational it
8 % 49
λ> let x = 2/7 in (2*x - x^2) / 3 :: Rational
8 % 49</pre>
=={{header|Java}}==
This example displays p-adic numbers in standard mathematical format, consisting of a possibly infinite list of digits extending leftwards from the p-adic point. p-adic numbers are given correct to O(prime^40) and the rational reconstruction is correct to O(prime^20).
<syntaxhighlight lang="java">
import java.util.ArrayList;
import java.util.Collections;
import java.util.List;
import java.util.stream.Collectors;
public final class PAdicNumbersBasic {
public static void main(String[] args) {
System.out.println("3-adic numbers:");
Padic padicOne = new Padic(3, -5, 9);
System.out.println("-5 / 9 => " + padicOne);
Padic padicTwo = new Padic(3, 47, 12);
System.out.println("47 / 12 => " + padicTwo);
Padic sum = padicOne.add(padicTwo);
System.out.println("sum => " + sum);
System.out.println("Rational = " + sum.convertToRational());
System.out.println();
System.out.println("7-adic numbers:");
padicOne = new Padic(7, 5, 8);
System.out.println("5 / 8 => " + padicOne);
padicTwo = new Padic(7, 353, 30809);
System.out.println("353 / 30809 => " + padicTwo);
sum = padicOne.add(padicTwo);
System.out.println("sum => " + sum);
System.out.println("Rational = " + sum.convertToRational());
}
}
final class Padic {
/**
* Create a p-adic, with p = aPrime, from the given rational 'aNumerator' / 'aDenominator'.
*/
public Padic(int aPrime, int aNumerator, int aDenominator) {
if ( aDenominator == 0 ) {
throw new IllegalArgumentException("Denominator cannot be zero");
}
prime = aPrime;
digits = new ArrayList<Integer>(DIGITS_SIZE);
order = 0;
// Process rational zero
if ( aNumerator == 0 ) {
order = MAX_ORDER;
return;
}
// Remove multiples of 'prime' and adjust the order of the p-adic number accordingly
while ( Math.floorMod(aNumerator, prime) == 0 ) {
aNumerator /= prime;
order += 1;
}
while ( Math.floorMod(aDenominator, prime) == 0 ) {
aDenominator /= prime;
order -= 1;
}
// Standard calculation of p-adic digits
final long inverse = moduloInverse(aDenominator);
while ( digits.size() < DIGITS_SIZE ) {
final int digit = Math.floorMod(aNumerator * inverse, prime);
digits.addLast(digit);
aNumerator -= digit * aDenominator;
if ( aNumerator != 0 ) {
// The denominator is not a power of a prime
int count = 0;
while ( Math.floorMod(aNumerator, prime) == 0 ) {
aNumerator /= prime;
count += 1;
}
for ( int i = count; i > 1; i-- ) {
digits.addLast(0);
}
}
}
}
/**
* Return the sum of this p-adic number and the given p-adic number.
*/
public Padic add(Padic aOther) {
if ( prime != aOther.prime ) {
throw new IllegalArgumentException("Cannot add p-adic's with different primes");
}
List<Integer> result = new ArrayList<Integer>();
// Adjust the digits so that the p-adic points are aligned
for ( int i = 0; i < -order + aOther.order; i++ ) {
aOther.digits.addFirst(0);
}
for ( int i = 0; i < -aOther.order + order; i++ ) {
digits.addFirst(0);
}
// Standard digit by digit addition
int carry = 0;
for ( int i = 0; i < Math.min(digits.size(), aOther.digits.size()); i++ ) {
final int sum = digits.get(i) + aOther.digits.get(i) + carry;
final int remainder = Math.floorMod(sum, prime);
carry = ( sum >= prime ) ? 1 : 0;
result.addLast(remainder);
}
// Reverse the changes made to the digits
for ( int i = 0; i < -order + aOther.order; i++ ) {
aOther.digits.removeFirst();
}
for ( int i = 0; i < -aOther.order + order; i++ ) {
digits.removeFirst();
}
return new Padic(prime, result, allZeroDigits(result) ? MAX_ORDER : Math.min(order, aOther.order));
}
/**
* Return the Rational representation of this p-adic number.
*/
public Rational convertToRational() {
List<Integer> numbers = new ArrayList<Integer>(digits);
// Zero
if ( numbers.isEmpty() || allZeroDigits(numbers) ) {
return new Rational(0, 1);
}
// Positive integer
if ( order >= 0 && endsWith(numbers, 0) ) {
for ( int i = 0; i < order; i++ ) {
numbers.addFirst(0);
}
return new Rational(convertToDecimal(numbers), 1);
}
// Negative integer
if ( order >= 0 && endsWith(numbers, prime - 1) ) {
negateList(numbers);
for ( int i = 0; i < order; i++ ) {
numbers.addFirst(0);
}
return new Rational(-convertToDecimal(numbers), 1);
}
// Rational
Padic sum = new Padic(prime, digits, order);
Padic self = new Padic(prime, digits, order);
int denominator = 1;
do {
sum = sum.add(self);
denominator += 1;
} while ( ! ( endsWith(sum.digits, 0) || endsWith(sum.digits, prime - 1) ) );
final boolean negative = endsWith(sum.digits, prime - 1);
if ( negative ) {
negateList(sum.digits);
}
int numerator = negative ? -convertToDecimal(sum.digits) : convertToDecimal(sum.digits);
if ( order > 0 ) {
numerator *= Math.pow(prime, order);
}
if ( order < 0 ) {
denominator *= Math.pow(prime, -order);
}
return new Rational(numerator, denominator);
}
/**
* Return a string representation of this p-adic.
*/
public String toString() {
List<Integer> numbers = new ArrayList<Integer>(digits);
padWithZeros(numbers);
Collections.reverse(numbers);
String numberString = numbers.stream().map(String::valueOf).collect(Collectors.joining());
StringBuilder builder = new StringBuilder(numberString);
if ( order >= 0 ) {
for ( int i = 0; i < order; i++ ) {
builder.append("0");
}
builder.append(".0");
} else {
builder.insert(builder.length() + order, ".");
while ( builder.toString().endsWith("0") ) {
builder.deleteCharAt(builder.length() - 1);
}
}
return " ..." + builder.toString().substring(builder.length() - PRECISION - 1);
}
// PRIVATE //
/**
* Create a p-adic, with p = 'aPrime', directly from a list of digits.
*
* With 'aOrder' = 0, the list [1, 2, 3, 4, 5] creates the p-adic ...54321.0
* 'aOrder' > 0 shifts the list 'aOrder' places to the left and
* 'aOrder' < 0 shifts the list 'aOrder' places to the right.
*/
private Padic(int aPrime, List<Integer> aDigits, int aOrder) {
prime = aPrime;
digits = new ArrayList<Integer>(aDigits);
order = aOrder;
}
/**
* Return the multiplicative inverse of the given decimal number modulo 'prime'.
*/
private int moduloInverse(int aNumber) {
int inverse = 1;
while ( Math.floorMod(inverse * aNumber, prime) != 1 ) {
inverse += 1;
}
return inverse;
}
/**
* Transform the given list of digits representing a p-adic number
* into a list which represents the negation of the p-adic number.
*/
private void negateList(List<Integer> aDigits) {
aDigits.set(0, Math.floorMod(prime - aDigits.get(0), prime));
for ( int i = 1; i < aDigits.size(); i++ ) {
aDigits.set(i, prime - 1 - aDigits.get(i));
}
}
/**
* Return the given list of base 'prime' integers converted to a decimal integer.
*/
private int convertToDecimal(List<Integer> aNumbers) {
int decimal = 0;
int multiple = 1;
for ( int number : aNumbers ) {
decimal += number * multiple;
multiple *= prime;
}
return decimal;
}
/**
* Return whether the given list consists of all zeros.
*/
private static boolean allZeroDigits(List<Integer> aList) {
return aList.stream().allMatch( i -> i == 0 );
}
/**
* The given list is padded on the right by zeros up to a maximum length of 'PRECISION'.
*/
private static void padWithZeros(List<Integer> aList) {
while ( aList.size() < DIGITS_SIZE ) {
aList.addLast(0);
}
}
/**
* Return whether the given list ends with multiple instances of the given number.
*/
private static boolean endsWith(List<Integer> aDigits, int aDigit) {
for ( int i = aDigits.size() - 1; i >= aDigits.size() - PRECISION / 2; i-- ) {
if ( aDigits.get(i) != aDigit ) {
return false;
}
}
return true;
}
private static class Rational {
public Rational(int aNumerator, int aDenominator) {
if ( aDenominator < 0 ) {
numerator = -aNumerator;
denominator = -aDenominator;
} else {
numerator = aNumerator;
denominator = aDenominator;
}
if ( aNumerator == 0 ) {
denominator = 1;
}
final int gcd = gcd(numerator, denominator);
numerator /= gcd;
denominator /= gcd;
}
public String toString() {
return numerator + " / " + denominator;
}
private int gcd(int aOne, int aTwo) {
if ( aTwo == 0 ) {
return Math.abs(aOne);
}
return gcd(aTwo, Math.floorMod(aOne, aTwo));
}
private int numerator;
private int denominator;
}
private List<Integer> digits;
private int order;
private final int prime;
private static final int MAX_ORDER = 1_000;
private static final int PRECISION = 40;
private static final int DIGITS_SIZE = PRECISION + 5;
}
</syntaxhighlight>
{{ out }}
<pre>
3-adic numbers:
-5 / 9 => ...22222222222222222222222222222222222222.11
47 / 12 => ...020202020202020202020202020202020202101.2
sum => ...20202020202020202020202020202020202101.01
Rational = 121 / 36
7-adic numbers:
5 / 8 => ...424242424242424242424242424242424242425.0
353 / 30809 => ...560462505550343461155520004023663643455.0
sum => ...315035233123101033613062431266421216213.0
Rational = 156869 / 246472
</pre>
=={{header|Julia}}==
Uses the Nemo abstract algebra library. The Nemo library's rational reconstruction function gives up quite easily,
so another alternative to FreeBasic's crat() using vector products is below.
<syntaxhighlight lang="julia">using Nemo, LinearAlgebra
set_printing_mode(FlintPadicField, :terse)
""" convert to Rational (rational reconstruction) """
function toRational(pa::padic)
rat = lift(QQ, pa)
r, den = BigInt(numerator(rat)), Int(denominator(rat))
p, k = Int(prime(parent(pa))), Int(precision(pa))
N = BigInt(p^k)
a1, a2 = [N, 0], [r, 1]
while dot(a1, a1) > dot(a2, a2)
q = dot(a1, a2) // dot(a2, a2)
a1, a2 = a2, a1 - BigInt(round(q)) * a2
end
if dot(a1, a1) < N
return (Rational{Int}(a1[1]) // Rational{Int}(a1[2])) // Int(den)
else
return Int(r) // den
end
end
function dstring(pa::padic)
u, v, n, p, k = pa.u, pa.v, pa.N, pa.parent.p, pa.parent.prec_max
d = digits(v > 0 ? u * p^v : u, base=pa.parent.p, pad=k)
return prod([i == k + v && v != 0 ? "$x . " : "$x " for (i, x) in enumerate(reverse(d))])
end
const DATA = [
[2, 1, 2, 4, 1, 1],
[4, 1, 2, 4, 3, 1],
[4, 1, 2, 5, 3, 1],
[4, 9, 5, 4, 8, 9],
[26, 25, 5, 4, -109, 125],
[49, 2, 7, 6, -4851, 2],
[-9, 5, 3, 8, 27, 7],
[5, 19, 2, 12, -101, 384],
# Base 10 10-adic p-adics are not allowed by Nemo library -- p must be a prime
# familiar digits
[11, 4, 2, 43, 679001, 207],
[-8, 9, 23, 9, 302113, 92],
[-22, 7, 3, 23, 46071, 379],
[-22, 7, 32749, 3, 46071, 379],
[35, 61, 5, 20, 9400, 109],
[-101, 109, 61, 7, 583376, 6649],
[-25, 26, 7, 13, 5571, 137],
[1, 4, 7, 11, 9263, 2837],
[122, 407, 7, 11, -517, 1477],
# more subtle
[5, 8, 7, 11, 353, 30809],
]
for (num1, den1, P, K, num2, den2) in DATA
Qp = PadicField(P, K)
a = Qp(QQ(num1 // den1))
b = Qp(QQ(num2 // den2))
c = a + b
r = toRational(c)
println(a, "\n", dstring(a), "\n", b, "\n", dstring(b), "\n+ =\n", c, "\n", dstring(c), " $r\n")
end
</syntaxhighlight>{{out}}
<pre>
2 + O(2^5)
0 0 1 0
1 + O(2^4)
0 0 0 1
+ =
3 + O(2^4)
0 0 1 1 3//1
4 + O(2^6)
0 1 0 0
3 + O(2^4)
0 0 1 1
+ =
7 + O(2^4)
0 1 1 1 -1//1
4 + O(2^7)
0 0 1 0 0
3 + O(2^5)
0 0 0 1 1
+ =
7 + O(2^5)
0 0 1 1 1 7//1
556 + O(5^4)
4 2 1 1
487 + O(5^4)
3 4 2 2
+ =
418 + O(5^4)
3 1 3 3 4//3
26/25 + O(5^2)
0 1 . 0 1
516/125 + O(5^1)
4 . 0 3 1
+ =
21/125 + O(5^1)
0 . 0 4 1 21//125
58849 + O(7^6)
3 3 3 4 0 0
56399 + O(7^6)
3 2 3 3 0 0
+ =
115248 + O(7^6)
6 6 0 0 0 0 0//1
5247 + O(3^8)
2 1 0 1 2 1 0 0
3753 + O(3^8)
1 2 0 1 1 0 0 0
+ =
2439 + O(3^8)
1 0 1 0 0 1 0 0 72//35
647 + O(2^12)
0 0 1 0 1 0 0 0 0 1 1 1
2697/128 + O(2^5)
1 0 1 0 1 . 0 0 0 1 0 0 1
+ =
3593/128 + O(2^5)
1 1 1 0 0 . 0 0 0 1 0 0 1 3593//128
11/4 + O(2^41)
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 . 1 1
2379619371607 + O(2^43)
0 1 0 0 0 1 0 1 0 1 0 0 0 0 0 1 1 0 0 0 1 0 1 1 1 1 0 0 0 0 0 1 0 1 0 0 1 0 1 0 1 1 1
+ =
722384464231/4 + O(2^41)
0 0 0 1 0 1 0 1 0 0 0 0 0 1 1 0 0 0 1 0 1 1 1 1 0 0 0 0 0 1 0 1 0 0 1 0 1 1 0 0 1 . 1 1 330311//1618052
200128073495 + O(23^9)
2 12 17 20 10 5 2 12 17
450288240894/23 + O(23^8)
5 17 5 17 6 0 10 12 . 2
+ =
1450928608353/23 + O(23^8)
18 12 3 4 11 3 0 6 . 2 1450928608353//23
40347076637 + O(3^23)
1 0 2 1 2 0 1 0 2 1 2 0 1 0 2 1 2 0 1 0 2 0 2
69303290076 + O(3^23)
2 0 1 2 1 2 1 2 2 1 2 1 0 0 2 2 0 1 1 2 1 0 0
+ =
15507187886 + O(3^23)
0 1 1 1 1 0 0 0 2 0 1 1 1 1 2 0 2 2 0 0 0 0 2 15507187886//1
30105603673496 + O(32749^3)
28070 18713 23389
4819014836161 + O(32749^3)
4493 8727 10145
+ =
34924618509657 + O(32749^3)
32563 27441 785 314159//2653
51592217117060 + O(5^20)
2 3 2 3 0 2 4 1 3 3 0 0 4 0 2 2 1 2 2 0
64744861847850 + O(5^20)
3 1 4 4 1 2 3 4 4 3 4 1 1 3 1 1 2 4 0 0
+ =
20969647324285 + O(5^20)
1 0 2 2 2 0 3 1 3 1 4 2 0 3 3 3 4 1 2 0 577215//6649
1701117681882 + O(61^7)
33 1 7 16 48 7 50
1701117687235/61 + O(61^6)
33 1 7 16 49 34 . 35
+ =
1758782693344/61 + O(61^6)
34 8 24 3 57 23 . 35 1758782693344//61
40991504402 + O(7^13)
2 6 5 0 5 4 4 0 1 6 1 2 2
46676457609 + O(7^13)
3 2 4 1 4 5 4 2 2 5 5 3 5
+ =
87667962011 + O(7^13)
6 2 2 2 3 3 1 2 4 4 6 6 0 141421//3562
494331686 + O(7^11)
1 5 1 5 1 5 1 5 1 5 2
1936205041 + O(7^11)
6 5 6 6 0 3 2 0 4 4 1
+ =
453209984 + O(7^11)
1 4 1 4 2 1 3 5 6 2 3 39889//11348
1778136580 + O(7^11)
6 2 0 3 0 6 2 4 4 4 3
384219886/7 + O(7^10)
1 2 3 4 3 5 4 6 4 1 . 1
+ =
967215488/7 + O(7^10)
3 2 6 5 3 1 2 4 1 4 . 1 967215488//7
1235829215 + O(7^11)
4 2 4 2 4 2 4 2 4 2 5
726006041 + O(7^11)
2 3 6 6 3 6 4 3 4 5 5
+ =
1961835256 + O(7^11)
6 6 4 2 1 2 1 6 2 1 3 -25145//36122
</pre>
=={{header|Nim}}==
{{trans|Go}}
Translation of Go with some modifications, especially using exceptions when an error is encountered.
<syntaxhighlight lang="nim">import math, strformat
const
Emx = 64 # Exponent maximum.
Dmx = 100000 # Approximation loop maximum.
Amx = 1048576 # Argument maximum.
PMax = 32749 # Prime maximum.
type
Ratio = tuple[a, b: int]
Padic = object
p: int # Prime.
k: int # Precision.
v: int
d: array[-Emx..(Emx-1), int]
PadicError = object of ValueError
proc r2pa(pa: var Padic; q: Ratio; sw: bool) =
## Convert "q" to p-adic number, set "sw" to print.
var (a, b) = q
if b == 0:
raise newException(PadicError, &"Wrong rational: {a}/{b}" )
if b < 0:
b = -b
a = -a
if abs(a) > Amx or b > Amx:
raise newException(PadicError, &"Rational exceeding limits: {a}/{b}")
if pa.p < 2:
raise newException(PadicError, &"Wrong value for p: {pa.p}")
if pa.k < 1:
raise newException(PadicError, &"Wrong value for k: {pa.k}")
pa.p = min(pa.p, PMax) # Maximum short prime.
pa.k = min(pa.k, Emx - 1) # Maximum array length.
if sw: echo &"{a}/{b} + 0({pa.p}^{pa.k})"
# Initialize.
pa.v = 0
pa.d.reset()
if a == 0: return
var i = 0
# Find -exponent of "p" in "b".
while b mod pa.p == 0:
b = b div pa.p
dec i
var s = 0
var r = b mod pa.p
# Modular inverse for small "p".
var b1 = 1
while b1 < pa.p:
inc s, r
if s >= pa.p: dec s, pa.p
if s == 1: break
inc b1
if b1 == pa.p:
raise newException(PadicError, "Impossible to compute inverse modulo")
pa.v = Emx
while true:
# Find exponent of "p" in "a".
while a mod pa.p == 0:
a = a div pa.p
inc i
# Valuation.
if pa.v == Emx: pa.v = i
# Upper bound.
if i >= Emx: break
# Check precision.
if i - pa.v > pa.k: break
# Next digit.
pa.d[i] = floorMod(a * b1, pa.p)
# Remainder - digit * divisor.
dec a, pa.d[i] * b
if a == 0: break
func dsum(pa: Padic): int =
## Horner's rule.
let t = min(pa.v, 0)
for i in countdown(pa.k - 1 + t, t):
var r = result
result *= pa.p
if r != 0 and (result div r - pa.p) != 0:
return -1 # Overflow.
inc result, pa.d[i]
func `+`(pa, pb: Padic): Padic =
## Add two p-adic numbers.
assert pa.p == pb.p and pa.k == pb.k
result.p = pa.p
result.k = pa.k
var c = 0
result.v = min(pa.v, pb.v)
for i in result.v..(pa.k + result.v):
inc c, pa.d[i] + pb.d[i]
if c >= pa.p:
result.d[i] = c - pa.p
c = 1
else:
result.d[i] = c
c = 0
func cmpt(pa: Padic): Padic =
## Return the complement.
var c = 1
result.p = pa.p
result.k = pa.k
result.v = pa.v
for i in pa.v..(pa.k + pa.v):
inc c, pa.p - 1 - pa.d[i]
if c >= pa.p:
result.d[i] = c - pa.p
c = 1
else:
result.d[i] = c
c = 0
func crat(pa: Padic): string =
## Rational reconstruction.
var s = pa
# Denominator count.
var i = 1
var fl = false
while i <= Dmx:
# Check for integer.
var j = pa.k - 1 + pa.v
while j >= pa.v:
if s.d[j] != 0: break
dec j
fl = (j - pa.v) * 2 < pa.k
if fl:
fl = false
break
# Check negative integer.
j = pa.k - 1 + pa.v
while j >= pa.v:
if pa.p - 1 - s.d[j] != 0: break
dec j
fl = (j - pa.v) * 2 < pa.k
if fl: break
# Repeatedly add "pa" to "s".
s = s + pa
inc i
if fl: s = s.cmpt()
# Numerator: weighted digit sum.
var x = s.dsum()
var y = i
if x < 0 or y > Dmx:
raise newException(PadicError, &"Error during rational reconstruction: {x}, {y}")
# Negative powers.
for i in pa.v..(-1): y *= pa.p
# Negative rational.
if fl: x = -x
result = $x
if y > 1: result.add &"/{y}"
func `$`(pa: Padic): string =
## String representation.
let t = min(pa.v, 0)
for i in countdown(pa.k - 1 + t, t):
result.add $pa.d[i]
if i == 0 and pa.v < 0: result.add "."
result.add " "
proc print(pa: Padic; sw: int) =
echo pa
# Rational approximation.
if sw != 0: echo pa.crat()
when isMainModule:
# Rational reconstruction depends on the precision
# until the dsum-loop overflows.
const Data = [[2, 1, 2, 4, 1, 1],
[4, 1, 2, 4, 3, 1],
[4, 1, 2, 5, 3, 1],
[4, 9, 5, 4, 8, 9],
[26, 25, 5, 4, -109, 125],
[49, 2, 7, 6, -4851, 2],
[-9, 5, 3, 8, 27, 7],
[5, 19, 2, 12, -101, 384],
# Two decadic pairs.
[2, 7, 10, 7, -1, 7],
[34, 21, 10, 9, -39034, 791],
# Familiar digits.
[11, 4, 2, 43, 679001, 207],
[-8, 9, 23, 9, 302113, 92],
[-22, 7, 3, 23, 46071, 379],
[-22, 7, 32749, 3, 46071, 379],
[35, 61, 5, 20, 9400, 109],
[-101, 109, 61, 7, 583376, 6649],
[-25, 26, 7, 13, 5571, 137],
[1, 4, 7, 11, 9263, 2837],
[122, 407, 7, 11, -517, 1477],
# More subtle.
[5, 8, 7, 11, 353, 30809]]
for d in Data:
try:
var a, b = Padic(p: d[2], k: d[3])
r2pa(a, (d[0], d[1]), true)
print(a, 0)
r2pa(b, (d[4], d[5]), true)
print(b, 0)
echo "+ ="
print(a + b, 1)
echo ""
except PadicError:
echo getCurrentExceptionMsg()</syntaxhighlight>
{{out}}
<pre>2/1 + 0(2^4)
0 0 1 0
1/1 + 0(2^4)
Line 1,096 ⟶ 2,647:
=={{header|Phix}}==
{{libheader|Phix/Class}}
<
constant EMX = 64 // exponent maximum (if indexing starts at -EMX)
constant DMX = 1e5 // approximation loop maximum
Line 1,363 ⟶ 2,914:
end if
printf(1,"\n")
end for</
{{out}}
<pre>
Line 1,385 ⟶ 2,936:
353/30809 + O(7^11) 2 3 6 6 3 6 4 3 4 5 5
47099/10977 + = 6 6 4 2 1 2 1 6 2 1 3
</pre>
=={{header|Raku}}==
<syntaxhighlight lang="raku" line># 20210225 Raku programming solution
#!/usr/bin/env raku
class Padic { has ($.p is default(2), %.v is default({})) is rw ;
method r2pa (Rat $x is copy, \p, \d) { # Reference: math.stackexchange.com/a/1187037
self.p = p ;
$x += p**d if $x < 0 ; # complement
my $lowerest = 0;
my ($num,$den) = $x.nude;
while ($den % p) == 0 { $den /= p and $lowerest-- }
$x = $num / $den;
while +self.v < d {
my %d = ^p Z=> (( $x «-« ^p ) »/» p )».&{ .denominator % p }; # .kv
for %d.keys { self.v.{$lowerest++} = $_ and last if %d{$_} != 0 }
$x = ($x - self.v.{$lowerest-1}) / p ;
}
self
}
method add (Padic \x, \d) {
my $div = 0;
my $lowerest = (self.v.keys.sort({.Int}).first,
x.v.keys.sort({.Int}).first ).min ;
return Padic.new:
p => self.p,
v => gather for ^d {
my $power = $lowerest + $_;
given ((self.v.{$power}//0)+(x.v.{$power}//0)+$div).polymod(x.p)
{ take ($power, .[0]).Slip and $div = .[1] }
}
}
method gist {
# en.wikipedia.org/wiki/P-adic_number#Notation
# my %H = (0..9) Z=> ('₀'..'₉'); # (0x2080 .. 0x2089);
# '⋯ ' ~ self.v ~ ' ' ~ [~] self.p.comb».&{ %H{$_} }
# express as a series
my %H = ( 0…9 ,'-') Z=> ( '⁰','¹','²','³','⁴'…'⁹','⁻');
[~] self.v.keys.sort({.Int}).map: {
' + ' ~ self.v.{$_} ~ '*' ~ self.p ~ [~] $_.comb».&{ %H{$_}} }
}
}
my @T;
for my \D = (
#`[[ these are not working
< 26/25 -109/125 5 4 >,
< 6/7 -5/7 10 7 >,
< 2/7 -3/7 10 7 >,
< 2/7 -1/7 10 7 >,
< 34/21 -39034/791 10 9 >,
#]]
#`[[[[[ Works
< 11/4 679001/207 2 43>,
< 11/4 679001/207 3 27 >,
< 5/19 -101/384 2 12>,
< -22/7 46071/379 7 13 >,
< -7/5 99/70 7 4> ,
< -101/109 583376/6649 61 7>,
< 122/407 -517/1477 7 11>,
< 2/1 1/1 2 4>,
< 4/1 3/1 2 4>,
< 4/1 3/1 2 5>,
< 4/9 8/9 5 4>,
< 11/4 679001/207 11 13 >,
< 1/4 9263/2837 7 11 >,
< 49/2 -4851/2 7 6 >,
< -9/5 27/7 3 8>,
< -22/7 46071/379 2 37 >,
< -22/7 46071/379 3 23 >,
< -101/109 583376/6649 2 40>,
< -101/109 583376/6649 32749 3>,
< -25/26 5571/137 7 13>,
#]]]]]
< 5/8 353/30809 7 11 >,
) -> \D {
given @T[0] = Padic.new { say D[0]~' = ', .r2pa: D[0],D[2],D[3] }
given @T[1] = Padic.new { say D[1]~' = ', .r2pa: D[1],D[2],D[3] }
given @T[0].add: @T[1], D[3] {
given @T[2] = Padic.new { .r2pa: D[0]+D[1], D[2], D[3] }
say "Addition result = ", $_.gist; #
unless ( $_.v.Str eq @T[2].v.Str ) {
say 'but ' ~ (D[0]+D[1]).nude.join('/') ~ ' = ' ~ @T[2].gist
}
}
}
</syntaxhighlight>
{{out}}
<pre>
5/8 = + 5*7⁰ + 2*7¹ + 4*7² + 2*7³ + 4*7⁴ + 2*7⁵ + 4*7⁶ + 2*7⁷ + 4*7⁸ + 2*7⁹ + 4*7¹⁰
353/30809 = + 5*7⁰ + 5*7¹ + 4*7² + 3*7³ + 4*7⁴ + 6*7⁵ + 3*7⁶ + 6*7⁷ + 6*7⁸ + 3*7⁹ + 2*7¹⁰
Addition result = + 3*7⁰ + 1*7¹ + 2*7² + 6*7³ + 1*7⁴ + 2*7⁵ + 1*7⁶ + 2*7⁷ + 4*7⁸ + 6*7⁹ + 6*7¹⁰
</pre>
Line 1,390 ⟶ 3,044:
{{trans|FreeBASIC}}
{{libheader|Wren-dynamic}}
<syntaxhighlight lang="wren">import "./dynamic" for Struct
// constants
Line 1,430 ⟶ 3,082:
if (a.abs > AMX || b > AMX) return -1
if (P < 2 || K < 1) return 1
P =
K =
if (sw != 0) {
System.write("%(a)/%(b) + ") // numerator, denominator
Line 1,493 ⟶ 3,145:
// Horner's rule
dsum() {
var t =
var s = 0
for (i in K - 1 + t..t) {
Line 1,567 ⟶ 3,219:
// print expansion
printf(sw) {
var t =
for (i in K - 1 + t..t) {
System.write(_d[i + EMX])
Line 1,582 ⟶ 3,234:
var c = 0
var r = Padic.new()
r.v =
for (i in r.v..K + r.v) {
c = c + _d[i+EMX] + b.d[i+EMX]
Line 1,665 ⟶ 3,317:
}
System.print()
}</
{{out}}
| |||