Mathematica code for computing the $p$-adic expansion of rational numbers

Question

Does anyone know any Mathematica code for computing the $p$-adic expansion of rational numbers? I.e. given a rational number $a/b,~a,b\in \mathbb{Z}$ and a prime number $p$, then compute the $p$-adic expansion of $a/b$ (to some order).

Answer 1

Something like this code should work:

pExpand[x : (_Rational | _Integer), p_?PrimeQ, n_ /; Positive[n]] := 
  Module[{q = p^n, num, den, v, e}, v = IntegerExponent[#, p] &
  /@ ({num, den} = #[x] & /@ {Numerator, Denominator});
  If[(e = v[[1]] - v[[2]]) > 0, num /= p^v[[1]], den /= p^v[[2]]];
  {e, Reverse@IntegerDigits[Mod[num PowerMod[den, -1, q], q], p, n]}];

The last line gives the $p$-adic valuation and digits. For example:

  pExpand[21/5, 7, 6] == {1, {2, 4, 5, 2, 1, 4}}

means that $\, 21/5 = 7^1\cdot (2 + 4\cdot7 + 5\cdot7^2 + 2\cdot7^3 + 1\cdot7^4 + 4\cdot7^5 + O(7^6)).$