Arithmetic derivative: Difference between revisions
(Initial Haskell version.) |
Thundergnat (talk | contribs) m (syntax highlighting fixup automation) |
||
Line 39: | Line 39: | ||
=={{header|C++}}== |
=={{header|C++}}== |
||
{{libheader|Boost}} |
{{libheader|Boost}} |
||
< |
<syntaxhighlight lang=cpp>#include <iomanip> |
||
#include <iostream> |
#include <iostream> |
||
Line 89: | Line 89: | ||
<< ") / 7 = " << arithmetic_derivative(p) / 7 << '\n'; |
<< ") / 7 = " << arithmetic_derivative(p) / 7 << '\n'; |
||
} |
} |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 140: | Line 140: | ||
Using ''float64'' (finessed a little) to avoid the unpleasantness of ''math/big'' for the stretch goal. |
Using ''float64'' (finessed a little) to avoid the unpleasantness of ''math/big'' for the stretch goal. |
||
Assumes that ''int'' type is 64 bit. |
Assumes that ''int'' type is 64 bit. |
||
< |
<syntaxhighlight lang=go>package main |
||
import ( |
import ( |
||
Line 185: | Line 185: | ||
fmt.Printf("D(10^%-2d) / 7 = %.0f\n", m, D(pow)/7) |
fmt.Printf("D(10^%-2d) / 7 = %.0f\n", m, D(pow)/7) |
||
} |
} |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 233: | Line 233: | ||
=={{header|Haskell}}== |
=={{header|Haskell}}== |
||
< |
<syntaxhighlight lang=haskell>import Control.Monad (forM_) |
||
import Data.List (intercalate) |
import Data.List (intercalate) |
||
import Data.List.Split (chunksOf) |
import Data.List.Split (chunksOf) |
||
Line 264: | Line 264: | ||
let q = 7 |
let q = 7 |
||
n = arithderiv (10^i) `quot` q |
n = arithderiv (10^i) `quot` q |
||
in printf "D(10^%-2d) / %d = %d\n" i q n</ |
in printf "D(10^%-2d) / %d = %d\n" i q n</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 313: | Line 313: | ||
=={{header|J}}== |
=={{header|J}}== |
||
Implementation: |
Implementation: |
||
< |
<syntaxhighlight lang=J>D=: {{ +/y%q:1>.|y }}"0</syntaxhighlight> |
||
In other words: find the sum of the argument divided by each of the sequence of prime factors of its absolute value (with a special case for zero -- we use the maximum of either 1 or that absolute value when finding the sequence of prime factors). |
In other words: find the sum of the argument divided by each of the sequence of prime factors of its absolute value (with a special case for zero -- we use the maximum of either 1 or that absolute value when finding the sequence of prime factors). |
||
Line 319: | Line 319: | ||
Task example: |
Task example: |
||
< |
<syntaxhighlight lang=J> D _99+i.20 10 |
||
_75 _77 _1 _272 _24 _49 _34 _96 _20 _123 |
_75 _77 _1 _272 _24 _49 _34 _96 _20 _123 |
||
_1 _140 _32 _45 _22 _124 _1 _43 _108 _176 |
_1 _140 _32 _45 _22 _124 _1 _43 _108 _176 |
||
Line 339: | Line 339: | ||
1 156 1 39 55 80 18 71 1 176 |
1 156 1 39 55 80 18 71 1 176 |
||
108 43 1 124 22 45 32 140 1 123 |
108 43 1 124 22 45 32 140 1 123 |
||
20 96 34 49 24 272 1 77 75 140</ |
20 96 34 49 24 272 1 77 75 140</syntaxhighlight> |
||
Also, it seems like it's worth verifying that order of evaluation does not create an ambiguity for the value of D: |
Also, it seems like it's worth verifying that order of evaluation does not create an ambiguity for the value of D: |
||
< |
<syntaxhighlight lang=J> 15 10 6 + 2 3 5 * D 15 10 6 |
||
31 31 31</ |
31 31 31</syntaxhighlight> |
||
Stretch task: |
Stretch task: |
||
< |
<syntaxhighlight lang=j> (D 10x^1+i.4 5)%7 |
||
1 20 300 4000 50000 |
1 20 300 4000 50000 |
||
600000 7000000 80000000 900000000 10000000000 |
600000 7000000 80000000 900000000 10000000000 |
||
110000000000 1200000000000 13000000000000 140000000000000 1500000000000000 |
110000000000 1200000000000 13000000000000 140000000000000 1500000000000000 |
||
16000000000000000 170000000000000000 1800000000000000000 19000000000000000000 200000000000000000000</ |
16000000000000000 170000000000000000 1800000000000000000 19000000000000000000 200000000000000000000</syntaxhighlight> |
||
=={{header|Julia}}== |
=={{header|Julia}}== |
||
< |
<syntaxhighlight lang=ruby>using Primes |
||
D(n) = n < 0 ? -D(-n) : n < 2 ? zero(n) : isprime(n) ? one(n) : typeof(n)(sum(e * n ÷ p for (p, e) in eachfactor(n))) |
D(n) = n < 0 ? -D(-n) : n < 2 ? zero(n) : isprime(n) ? one(n) : typeof(n)(sum(e * n ÷ p for (p, e) in eachfactor(n))) |
||
Line 365: | Line 365: | ||
println("D for 10^", rpad(m, 3), "divided by 7 is ", D(Int128(10)^m) ÷ 7) |
println("D for 10^", rpad(m, 3), "divided by 7 is ", D(Int128(10)^m) ÷ 7) |
||
end |
end |
||
</ |
</syntaxhighlight>{{out}} |
||
<pre> |
<pre> |
||
-75 -77 -1 -272 -24 -49 -34 -96 -20 -123 |
-75 -77 -1 -272 -24 -49 -34 -96 -20 -123 |
||
Line 413: | Line 413: | ||
{{trans|J}} |
{{trans|J}} |
||
{{libheader|ntheory}} |
{{libheader|ntheory}} |
||
< |
<syntaxhighlight lang=perl>use v5.36; |
||
use bigint; |
use bigint; |
||
no warnings 'uninitialized'; |
no warnings 'uninitialized'; |
||
Line 430: | Line 430: | ||
say table 10, map { D $_ } -99 .. 100; |
say table 10, map { D $_ } -99 .. 100; |
||
say join "\n", map { sprintf('D(10**%-2d) / 7 == ', $_) . D(10**$_) / 7 } 1 .. 20; |
say join "\n", map { sprintf('D(10**%-2d) / 7 == ', $_) . D(10**$_) / 7 } 1 .. 20; |
||
</syntaxhighlight> |
|||
</lang> |
|||
{{out}} |
{{out}} |
||
<pre> -75 -77 -1 -272 -24 -49 -34 -96 -20 -123 |
<pre> -75 -77 -1 -272 -24 -49 -34 -96 -20 -123 |
||
Line 475: | Line 475: | ||
=={{header|Phix}}== |
=={{header|Phix}}== |
||
<!--< |
<!--<syntaxhighlight lang=Phix>(phixonline)--> |
||
<span style="color: #008080;">with</span> <span style="color: #008080;">javascript_semantics</span> |
<span style="color: #008080;">with</span> <span style="color: #008080;">javascript_semantics</span> |
||
<span style="color: #008080;">include</span> <span style="color: #004080;">mpfr</span><span style="color: #0000FF;">.</span><span style="color: #000000;">e</span> |
<span style="color: #008080;">include</span> <span style="color: #004080;">mpfr</span><span style="color: #0000FF;">.</span><span style="color: #000000;">e</span> |
||
Line 516: | Line 516: | ||
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"D(10^%d)/7 = %s\n"</span><span style="color: #0000FF;">,{</span><span style="color: #000000;">m</span><span style="color: #0000FF;">,</span><span style="color: #7060A8;">mpz_get_str</span><span style="color: #0000FF;">(</span><span style="color: #000000;">n</span><span style="color: #0000FF;">)})</span> |
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"D(10^%d)/7 = %s\n"</span><span style="color: #0000FF;">,{</span><span style="color: #000000;">m</span><span style="color: #0000FF;">,</span><span style="color: #7060A8;">mpz_get_str</span><span style="color: #0000FF;">(</span><span style="color: #000000;">n</span><span style="color: #0000FF;">)})</span> |
||
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span> |
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span> |
||
<!--</ |
<!--</syntaxhighlight>--> |
||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 563: | Line 563: | ||
=={{header|Python}}== |
=={{header|Python}}== |
||
< |
<syntaxhighlight lang=python>from sympy.ntheory import factorint |
||
def D(n): |
def D(n): |
||
Line 583: | Line 583: | ||
print('(D for 10**{}) divided by 7 is {}'.format(m, D(10 ** m) // 7)) |
print('(D for 10**{}) divided by 7 is {}'.format(m, D(10 ** m) // 7)) |
||
</ |
</syntaxhighlight>{{out}} |
||
<pre> |
<pre> |
||
-75 -77 -1 -272 -24 -49 -34 -96 -20 -123 |
-75 -77 -1 -272 -24 -49 -34 -96 -20 -123 |
||
Line 629: | Line 629: | ||
=={{header|Raku}}== |
=={{header|Raku}}== |
||
< |
<syntaxhighlight lang=perl6>use Prime::Factor; |
||
multi D (0) { 0 } |
multi D (0) { 0 } |
||
Line 642: | Line 642: | ||
put ''; |
put ''; |
||
put join "\n", (1..20).map: { sprintf "D(10**%-2d) / 7 == %d", $_, D(10**$_) / 7 }</ |
put join "\n", (1..20).map: { sprintf "D(10**%-2d) / 7 == %d", $_, D(10**$_) / 7 }</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> -75 -77 -1 -272 -24 -49 -34 -96 -20 -123 |
<pre> -75 -77 -1 -272 -24 -49 -34 -96 -20 -123 |
||
Line 690: | Line 690: | ||
{{libheader|Wren-fmt}} |
{{libheader|Wren-fmt}} |
||
As integer arithmetic in Wren is inaccurate above 2^53 we need to use BigInt here. |
As integer arithmetic in Wren is inaccurate above 2^53 we need to use BigInt here. |
||
< |
<syntaxhighlight lang=ecmascript>import "./big" for BigInt |
||
import "./fmt" for Fmt |
import "./fmt" for Fmt |
||
Line 710: | Line 710: | ||
for (m in 1..20) { |
for (m in 1..20) { |
||
Fmt.print("D(10^$-2d) / 7 = $i", m, D.call(BigInt.ten.pow(m))/7) |
Fmt.print("D(10^$-2d) / 7 = $i", m, D.call(BigInt.ten.pow(m))/7) |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
Revision as of 01:11, 26 August 2022
The arithmetic derivative of an integer (more specifically, the Lagarias arithmetic derivative) is a function defined for integers, based on prime factorization, by analogy with the product rule for the derivative of a function that is used in mathematical analysis. Accordingly, for natural numbers n, the arithmetic derivative D(n) is defined as follows:
- D(0) = D(1) = 0.
- D(p) = 1 for any prime p.
- D(mn) = D(m)n + mD(n) for any m,n ∈ N. (Leibniz rule for derivatives).
Additionally, for negative integers the arithmetic derivative may be defined as -D(-n) (n < 0).
- Examples
D(2) = 1 and D(3) = 1 (both are prime) so if mn = 2 * 3, D(6) = (1)(3) + (1)(2) = 5.
D(9) = D(3)(3) + D(3)(3) = 6
D(27) = D(3)*9 + D(9)*3 = 9 + 18 = 27
D(30) = D(5)(6) + D(6)(5) = 6 + 5 * 5 = 31.
- Task
Find and show the arithmetic derivatives for -99 through 100.
- Stretch task
Find (the arithmetic derivative of 10^m) then divided by 7, where m is from 1 to 20.
- See also
C++
#include <iomanip>
#include <iostream>
#include <boost/multiprecision/cpp_int.hpp>
template <typename IntegerType>
IntegerType arithmetic_derivative(IntegerType n) {
bool negative = n < 0;
if (negative)
n = -n;
if (n < 2)
return 0;
IntegerType sum = 0, count = 0, m = n;
while ((m & 1) == 0) {
m >>= 1;
count += n;
}
if (count > 0)
sum += count / 2;
for (IntegerType p = 3, sq = 9; sq <= m; p += 2) {
count = 0;
while (m % p == 0) {
m /= p;
count += n;
}
if (count > 0)
sum += count / p;
sq += (p + 1) << 2;
}
if (m > 1)
sum += n / m;
if (negative)
sum = -sum;
return sum;
}
int main() {
using boost::multiprecision::int128_t;
for (int n = -99, i = 0; n <= 100; ++n, ++i) {
std::cout << std::setw(4) << arithmetic_derivative(n)
<< ((i + 1) % 10 == 0 ? '\n' : ' ');
}
int128_t p = 10;
std::cout << '\n';
for (int i = 0; i < 20; ++i, p *= 10) {
std::cout << "D(10^" << std::setw(2) << i + 1
<< ") / 7 = " << arithmetic_derivative(p) / 7 << '\n';
}
}
- Output:
-75 -77 -1 -272 -24 -49 -34 -96 -20 -123 -1 -140 -32 -45 -22 -124 -1 -43 -108 -176 -1 -71 -18 -80 -55 -39 -1 -156 -1 -59 -26 -72 -1 -61 -18 -192 -51 -33 -1 -92 -1 -31 -22 -92 -16 -81 -1 -56 -20 -45 -14 -112 -1 -25 -39 -48 -1 -41 -1 -68 -16 -21 -1 -60 -12 -19 -14 -80 -1 -31 -1 -32 -27 -15 -10 -44 -1 -13 -10 -24 -1 -21 -1 -32 -8 -9 -1 -16 -1 -7 -6 -12 -1 -5 -1 -4 -1 -1 0 0 0 1 1 4 1 5 1 12 6 7 1 16 1 9 8 32 1 21 1 24 10 13 1 44 10 15 27 32 1 31 1 80 14 19 12 60 1 21 16 68 1 41 1 48 39 25 1 112 14 45 20 56 1 81 16 92 22 31 1 92 1 33 51 192 18 61 1 72 26 59 1 156 1 39 55 80 18 71 1 176 108 43 1 124 22 45 32 140 1 123 20 96 34 49 24 272 1 77 75 140 D(10^ 1) / 7 = 1 D(10^ 2) / 7 = 20 D(10^ 3) / 7 = 300 D(10^ 4) / 7 = 4000 D(10^ 5) / 7 = 50000 D(10^ 6) / 7 = 600000 D(10^ 7) / 7 = 7000000 D(10^ 8) / 7 = 80000000 D(10^ 9) / 7 = 900000000 D(10^10) / 7 = 10000000000 D(10^11) / 7 = 110000000000 D(10^12) / 7 = 1200000000000 D(10^13) / 7 = 13000000000000 D(10^14) / 7 = 140000000000000 D(10^15) / 7 = 1500000000000000 D(10^16) / 7 = 16000000000000000 D(10^17) / 7 = 170000000000000000 D(10^18) / 7 = 1800000000000000000 D(10^19) / 7 = 19000000000000000000 D(10^20) / 7 = 200000000000000000000
Go
Using float64 (finessed a little) to avoid the unpleasantness of math/big for the stretch goal. Assumes that int type is 64 bit.
package main
import (
"fmt"
"rcu"
)
func D(n float64) float64 {
if n < 0 {
return -D(-n)
}
if n < 2 {
return 0
}
var f []int
if n < 1e19 {
f = rcu.PrimeFactors(int(n))
} else {
g := int(n / 100)
f = rcu.PrimeFactors(g)
f = append(f, []int{2, 2, 5, 5}...)
}
c := len(f)
if c == 1 {
return 1
}
if c == 2 {
return float64(f[0] + f[1])
}
d := n / float64(f[0])
return D(d)*float64(f[0]) + d
}
func main() {
ad := make([]int, 200)
for n := -99; n < 101; n++ {
ad[n+99] = int(D(float64(n)))
}
rcu.PrintTable(ad, 10, 4, false)
fmt.Println()
pow := 1.0
for m := 1; m < 21; m++ {
pow *= 10
fmt.Printf("D(10^%-2d) / 7 = %.0f\n", m, D(pow)/7)
}
}
- Output:
-75 -77 -1 -272 -24 -49 -34 -96 -20 -123 -1 -140 -32 -45 -22 -124 -1 -43 -108 -176 -1 -71 -18 -80 -55 -39 -1 -156 -1 -59 -26 -72 -1 -61 -18 -192 -51 -33 -1 -92 -1 -31 -22 -92 -16 -81 -1 -56 -20 -45 -14 -112 -1 -25 -39 -48 -1 -41 -1 -68 -16 -21 -1 -60 -12 -19 -14 -80 -1 -31 -1 -32 -27 -15 -10 -44 -1 -13 -10 -24 -1 -21 -1 -32 -8 -9 -1 -16 -1 -7 -6 -12 -1 -5 -1 -4 -1 -1 0 0 0 1 1 4 1 5 1 12 6 7 1 16 1 9 8 32 1 21 1 24 10 13 1 44 10 15 27 32 1 31 1 80 14 19 12 60 1 21 16 68 1 41 1 48 39 25 1 112 14 45 20 56 1 81 16 92 22 31 1 92 1 33 51 192 18 61 1 72 26 59 1 156 1 39 55 80 18 71 1 176 108 43 1 124 22 45 32 140 1 123 20 96 34 49 24 272 1 77 75 140 D(10^1 ) / 7 = 1 D(10^2 ) / 7 = 20 D(10^3 ) / 7 = 300 D(10^4 ) / 7 = 4000 D(10^5 ) / 7 = 50000 D(10^6 ) / 7 = 600000 D(10^7 ) / 7 = 7000000 D(10^8 ) / 7 = 80000000 D(10^9 ) / 7 = 900000000 D(10^10) / 7 = 10000000000 D(10^11) / 7 = 110000000000 D(10^12) / 7 = 1200000000000 D(10^13) / 7 = 13000000000000 D(10^14) / 7 = 140000000000000 D(10^15) / 7 = 1500000000000000 D(10^16) / 7 = 16000000000000000 D(10^17) / 7 = 170000000000000000 D(10^18) / 7 = 1800000000000000000 D(10^19) / 7 = 19000000000000000000 D(10^20) / 7 = 200000000000000000000
Haskell
import Control.Monad (forM_)
import Data.List (intercalate)
import Data.List.Split (chunksOf)
import Math.NumberTheory.Primes (factorise, unPrime)
import Text.Printf (printf)
-- The arithmetic derivative of a number, which is assumed to be non-negative.
arithderiv_ :: Integer -> Integer
arithderiv_ 0 = 0
arithderiv_ n = foldr step 0 $ factorise n
where step (p, v) s = s + n `quot` unPrime p * fromIntegral v
-- The arithmetic derivative of any integer.
arithderiv :: Integer -> Integer
arithderiv n | n < 0 = negate $ arithderiv_ (negate n)
| otherwise = arithderiv_ n
printTable :: [Integer] -> IO ()
printTable = putStrLn
. intercalate "\n"
. map unwords
. chunksOf 10
. map (printf "%5d")
main :: IO ()
main = do
printTable [arithderiv n | n <- [-99..100]]
putStrLn ""
forM_ [1..20 :: Integer] $ \i ->
let q = 7
n = arithderiv (10^i) `quot` q
in printf "D(10^%-2d) / %d = %d\n" i q n
- Output:
$ arithderiv -75 -77 -1 -272 -24 -49 -34 -96 -20 -123 -1 -140 -32 -45 -22 -124 -1 -43 -108 -176 -1 -71 -18 -80 -55 -39 -1 -156 -1 -59 -26 -72 -1 -61 -18 -192 -51 -33 -1 -92 -1 -31 -22 -92 -16 -81 -1 -56 -20 -45 -14 -112 -1 -25 -39 -48 -1 -41 -1 -68 -16 -21 -1 -60 -12 -19 -14 -80 -1 -31 -1 -32 -27 -15 -10 -44 -1 -13 -10 -24 -1 -21 -1 -32 -8 -9 -1 -16 -1 -7 -6 -12 -1 -5 -1 -4 -1 -1 0 0 0 1 1 4 1 5 1 12 6 7 1 16 1 9 8 32 1 21 1 24 10 13 1 44 10 15 27 32 1 31 1 80 14 19 12 60 1 21 16 68 1 41 1 48 39 25 1 112 14 45 20 56 1 81 16 92 22 31 1 92 1 33 51 192 18 61 1 72 26 59 1 156 1 39 55 80 18 71 1 176 108 43 1 124 22 45 32 140 1 123 20 96 34 49 24 272 1 77 75 140 D(10^1 ) / 7 = 1 D(10^2 ) / 7 = 20 D(10^3 ) / 7 = 300 D(10^4 ) / 7 = 4000 D(10^5 ) / 7 = 50000 D(10^6 ) / 7 = 600000 D(10^7 ) / 7 = 7000000 D(10^8 ) / 7 = 80000000 D(10^9 ) / 7 = 900000000 D(10^10) / 7 = 10000000000 D(10^11) / 7 = 110000000000 D(10^12) / 7 = 1200000000000 D(10^13) / 7 = 13000000000000 D(10^14) / 7 = 140000000000000 D(10^15) / 7 = 1500000000000000 D(10^16) / 7 = 16000000000000000 D(10^17) / 7 = 170000000000000000 D(10^18) / 7 = 1800000000000000000 D(10^19) / 7 = 19000000000000000000 D(10^20) / 7 = 200000000000000000000
J
Implementation:
D=: {{ +/y%q:1>.|y }}"0
In other words: find the sum of the argument divided by each of the sequence of prime factors of its absolute value (with a special case for zero -- we use the maximum of either 1 or that absolute value when finding the sequence of prime factors).
Task example:
D _99+i.20 10
_75 _77 _1 _272 _24 _49 _34 _96 _20 _123
_1 _140 _32 _45 _22 _124 _1 _43 _108 _176
_1 _71 _18 _80 _55 _39 _1 _156 _1 _59
_26 _72 _1 _61 _18 _192 _51 _33 _1 _92
_1 _31 _22 _92 _16 _81 _1 _56 _20 _45
_14 _112 _1 _25 _39 _48 _1 _41 _1 _68
_16 _21 _1 _60 _12 _19 _14 _80 _1 _31
_1 _32 _27 _15 _10 _44 _1 _13 _10 _24
_1 _21 _1 _32 _8 _9 _1 _16 _1 _7
_6 _12 _1 _5 _1 _4 _1 _1 0 0
0 1 1 4 1 5 1 12 6 7
1 16 1 9 8 32 1 21 1 24
10 13 1 44 10 15 27 32 1 31
1 80 14 19 12 60 1 21 16 68
1 41 1 48 39 25 1 112 14 45
20 56 1 81 16 92 22 31 1 92
1 33 51 192 18 61 1 72 26 59
1 156 1 39 55 80 18 71 1 176
108 43 1 124 22 45 32 140 1 123
20 96 34 49 24 272 1 77 75 140
Also, it seems like it's worth verifying that order of evaluation does not create an ambiguity for the value of D:
15 10 6 + 2 3 5 * D 15 10 6
31 31 31
Stretch task:
(D 10x^1+i.4 5)%7
1 20 300 4000 50000
600000 7000000 80000000 900000000 10000000000
110000000000 1200000000000 13000000000000 140000000000000 1500000000000000
16000000000000000 170000000000000000 1800000000000000000 19000000000000000000 200000000000000000000
Julia
using Primes
D(n) = n < 0 ? -D(-n) : n < 2 ? zero(n) : isprime(n) ? one(n) : typeof(n)(sum(e * n ÷ p for (p, e) in eachfactor(n)))
foreach(p -> print(lpad(p[2], 5), p[1] % 10 == 0 ? "\n" : ""), pairs(map(D, -99:100)))
println()
for m in 1:20
println("D for 10^", rpad(m, 3), "divided by 7 is ", D(Int128(10)^m) ÷ 7)
end
- Output:
-75 -77 -1 -272 -24 -49 -34 -96 -20 -123 -1 -140 -32 -45 -22 -124 -1 -43 -108 -176 -1 -71 -18 -80 -55 -39 -1 -156 -1 -59 -26 -72 -1 -61 -18 -192 -51 -33 -1 -92 -1 -31 -22 -92 -16 -81 -1 -56 -20 -45 -14 -112 -1 -25 -39 -48 -1 -41 -1 -68 -16 -21 -1 -60 -12 -19 -14 -80 -1 -31 -1 -32 -27 -15 -10 -44 -1 -13 -10 -24 -1 -21 -1 -32 -8 -9 -1 -16 -1 -7 -6 -12 -1 -5 -1 -4 -1 -1 0 0 0 1 1 4 1 5 1 12 6 7 1 16 1 9 8 32 1 21 1 24 10 13 1 44 10 15 27 32 1 31 1 80 14 19 12 60 1 21 16 68 1 41 1 48 39 25 1 112 14 45 20 56 1 81 16 92 22 31 1 92 1 33 51 192 18 61 1 72 26 59 1 156 1 39 55 80 18 71 1 176 108 43 1 124 22 45 32 140 1 123 20 96 34 49 24 272 1 77 75 140 D for 10^1 divided by 7 is 1 D for 10^2 divided by 7 is 20 D for 10^3 divided by 7 is 300 D for 10^4 divided by 7 is 4000 D for 10^5 divided by 7 is 50000 D for 10^6 divided by 7 is 600000 D for 10^7 divided by 7 is 7000000 D for 10^8 divided by 7 is 80000000 D for 10^9 divided by 7 is 900000000 D for 10^10 divided by 7 is 10000000000 D for 10^11 divided by 7 is 110000000000 D for 10^12 divided by 7 is 1200000000000 D for 10^13 divided by 7 is 13000000000000 D for 10^14 divided by 7 is 140000000000000 D for 10^15 divided by 7 is 1500000000000000 D for 10^16 divided by 7 is 16000000000000000 D for 10^17 divided by 7 is 170000000000000000 D for 10^18 divided by 7 is 1800000000000000000 D for 10^19 divided by 7 is 19000000000000000000 D for 10^20 divided by 7 is 200000000000000000000
Perl
use v5.36;
use bigint;
no warnings 'uninitialized';
use List::Util 'max';
use ntheory 'factor';
sub table ($c, @V) { my $t = $c * (my $w = 2 + length max @V); ( sprintf( ('%'.$w.'d')x@V, @V) ) =~ s/.{1,$t}\K/\n/gr }
sub D ($n) {
my(%f, $s);
$f{$_}++ for factor max 1, my $nabs = abs $n;
map { $s += $nabs * $f{$_} / $_ } keys %f;
$n > 0 ? $s : -$s;
}
say table 10, map { D $_ } -99 .. 100;
say join "\n", map { sprintf('D(10**%-2d) / 7 == ', $_) . D(10**$_) / 7 } 1 .. 20;
- Output:
-75 -77 -1 -272 -24 -49 -34 -96 -20 -123 -1 -140 -32 -45 -22 -124 -1 -43 -108 -176 -1 -71 -18 -80 -55 -39 -1 -156 -1 -59 -26 -72 -1 -61 -18 -192 -51 -33 -1 -92 -1 -31 -22 -92 -16 -81 -1 -56 -20 -45 -14 -112 -1 -25 -39 -48 -1 -41 -1 -68 -16 -21 -1 -60 -12 -19 -14 -80 -1 -31 -1 -32 -27 -15 -10 -44 -1 -13 -10 -24 -1 -21 -1 -32 -8 -9 -1 -16 -1 -7 -6 -12 -1 -5 -1 -4 -1 -1 0 0 0 1 1 4 1 5 1 12 6 7 1 16 1 9 8 32 1 21 1 24 10 13 1 44 10 15 27 32 1 31 1 80 14 19 12 60 1 21 16 68 1 41 1 48 39 25 1 112 14 45 20 56 1 81 16 92 22 31 1 92 1 33 51 192 18 61 1 72 26 59 1 156 1 39 55 80 18 71 1 176 108 43 1 124 22 45 32 140 1 123 20 96 34 49 24 272 1 77 75 140 D(10**1 ) / 7 == 1 D(10**2 ) / 7 == 20 D(10**3 ) / 7 == 300 D(10**4 ) / 7 == 4000 D(10**5 ) / 7 == 50000 D(10**6 ) / 7 == 600000 D(10**7 ) / 7 == 7000000 D(10**8 ) / 7 == 80000000 D(10**9 ) / 7 == 900000000 D(10**10) / 7 == 10000000000 D(10**11) / 7 == 110000000000 D(10**12) / 7 == 1200000000000 D(10**13) / 7 == 13000000000000 D(10**14) / 7 == 140000000000000 D(10**15) / 7 == 1500000000000000 D(10**16) / 7 == 16000000000000000 D(10**17) / 7 == 170000000000000000 D(10**18) / 7 == 1800000000000000000 D(10**19) / 7 == 19000000000000000000 D(10**20) / 7 == 200000000000000000000
Phix
with javascript_semantics include mpfr.e procedure D(mpz n) integer s = mpz_cmp_si(n,0) if s<0 then mpz_neg(n,n) end if if mpz_cmp_si(n,2)<0 then mpz_set_si(n,0) else sequence f = mpz_prime_factors(n) integer c = sum(vslice(f,2)), f1 = f[1][1] if c=1 then mpz_set_si(n,1) elsif c=2 then mpz_set_si(n,f1 + iff(length(f)=1?f1:f[2][1])) else assert(mpz_fdiv_q_ui(n,n,f1)=0) mpz d = mpz_init_set(n) D(n) mpz_mul_si(n,n,f1) mpz_add(n,n,d) end if if s<0 then mpz_neg(n,n) end if end if end procedure sequence res = repeat(0,200) mpz n = mpz_init() for i=-99 to 100 do mpz_set_si(n,i) D(n) res[i+100] = mpz_get_str(n) end for printf(1,"%s\n\n",{join_by(res,1,10," ",fmt:="%4s")}) for m=1 to 20 do mpz_ui_pow_ui(n,10,m) D(n) assert(mpz_fdiv_q_ui(n,n,7)=0) printf(1,"D(10^%d)/7 = %s\n",{m,mpz_get_str(n)}) end for
- Output:
-75 -77 -1 -272 -24 -49 -34 -96 -20 -123 -1 -140 -32 -45 -22 -124 -1 -43 -108 -176 -1 -71 -18 -80 -55 -39 -1 -156 -1 -59 -26 -72 -1 -61 -18 -192 -51 -33 -1 -92 -1 -31 -22 -92 -16 -81 -1 -56 -20 -45 -14 -112 -1 -25 -39 -48 -1 -41 -1 -68 -16 -21 -1 -60 -12 -19 -14 -80 -1 -31 -1 -32 -27 -15 -10 -44 -1 -13 -10 -24 -1 -21 -1 -32 -8 -9 -1 -16 -1 -7 -6 -12 -1 -5 -1 -4 -1 -1 0 0 0 1 1 4 1 5 1 12 6 7 1 16 1 9 8 32 1 21 1 24 10 13 1 44 10 15 27 32 1 31 1 80 14 19 12 60 1 21 16 68 1 41 1 48 39 25 1 112 14 45 20 56 1 81 16 92 22 31 1 92 1 33 51 192 18 61 1 72 26 59 1 156 1 39 55 80 18 71 1 176 108 43 1 124 22 45 32 140 1 123 20 96 34 49 24 272 1 77 75 140 D(10^1)/7 = 1 D(10^2)/7 = 20 D(10^3)/7 = 300 D(10^4)/7 = 4000 D(10^5)/7 = 50000 D(10^6)/7 = 600000 D(10^7)/7 = 7000000 D(10^8)/7 = 80000000 D(10^9)/7 = 900000000 D(10^10)/7 = 10000000000 D(10^11)/7 = 110000000000 D(10^12)/7 = 1200000000000 D(10^13)/7 = 13000000000000 D(10^14)/7 = 140000000000000 D(10^15)/7 = 1500000000000000 D(10^16)/7 = 16000000000000000 D(10^17)/7 = 170000000000000000 D(10^18)/7 = 1800000000000000000 D(10^19)/7 = 19000000000000000000 D(10^20)/7 = 200000000000000000000
Python
from sympy.ntheory import factorint
def D(n):
if n < 0:
return -D(-n)
elif n < 2:
return 0
else:
fdict = factorint(n)
if len(fdict) == 1 and 1 in fdict: # is prime
return 1
return sum([n * e // p for p, e in fdict.items()])
for n in range(-99, 101):
print('{:5}'.format(D(n)), end='\n' if n % 10 == 0 else '')
print()
for m in range(1, 21):
print('(D for 10**{}) divided by 7 is {}'.format(m, D(10 ** m) // 7))
- Output:
-75 -77 -1 -272 -24 -49 -34 -96 -20 -123 -1 -140 -32 -45 -22 -124 -1 -43 -108 -176 -1 -71 -18 -80 -55 -39 -1 -156 -1 -59 -26 -72 -1 -61 -18 -192 -51 -33 -1 -92 -1 -31 -22 -92 -16 -81 -1 -56 -20 -45 -14 -112 -1 -25 -39 -48 -1 -41 -1 -68 -16 -21 -1 -60 -12 -19 -14 -80 -1 -31 -1 -32 -27 -15 -10 -44 -1 -13 -10 -24 -1 -21 -1 -32 -8 -9 -1 -16 -1 -7 -6 -12 -1 -5 -1 -4 -1 -1 0 0 0 1 1 4 1 5 1 12 6 7 1 16 1 9 8 32 1 21 1 24 10 13 1 44 10 15 27 32 1 31 1 80 14 19 12 60 1 21 16 68 1 41 1 48 39 25 1 112 14 45 20 56 1 81 16 92 22 31 1 92 1 33 51 192 18 61 1 72 26 59 1 156 1 39 55 80 18 71 1 176 108 43 1 124 22 45 32 140 1 123 20 96 34 49 24 272 1 77 75 140 (D for 10**1) divided by 7 is 1 (D for 10**2) divided by 7 is 20 (D for 10**3) divided by 7 is 300 (D for 10**4) divided by 7 is 4000 (D for 10**5) divided by 7 is 50000 (D for 10**6) divided by 7 is 600000 (D for 10**7) divided by 7 is 7000000 (D for 10**8) divided by 7 is 80000000 (D for 10**9) divided by 7 is 900000000 (D for 10**10) divided by 7 is 10000000000 (D for 10**11) divided by 7 is 110000000000 (D for 10**12) divided by 7 is 1200000000000 (D for 10**13) divided by 7 is 13000000000000 (D for 10**14) divided by 7 is 140000000000000 (D for 10**15) divided by 7 is 1500000000000000 (D for 10**16) divided by 7 is 16000000000000000 (D for 10**17) divided by 7 is 170000000000000000 (D for 10**18) divided by 7 is 1800000000000000000 (D for 10**19) divided by 7 is 19000000000000000000 (D for 10**20) divided by 7 is 200000000000000000000
Raku
use Prime::Factor;
multi D (0) { 0 }
multi D (1) { 0 }
multi D ($n where &is-prime) { 1 }
multi D ($n where * < 0 ) { -D -$n }
multi D ($n) { sum $n.&prime-factors.Bag.map: { $n × .value / .key } }
put (-99 .. 100).map(&D).batch(10)».fmt("%4d").join: "\n";
put '';
put join "\n", (1..20).map: { sprintf "D(10**%-2d) / 7 == %d", $_, D(10**$_) / 7 }
- Output:
-75 -77 -1 -272 -24 -49 -34 -96 -20 -123 -1 -140 -32 -45 -22 -124 -1 -43 -108 -176 -1 -71 -18 -80 -55 -39 -1 -156 -1 -59 -26 -72 -1 -61 -18 -192 -51 -33 -1 -92 -1 -31 -22 -92 -16 -81 -1 -56 -20 -45 -14 -112 -1 -25 -39 -48 -1 -41 -1 -68 -16 -21 -1 -60 -12 -19 -14 -80 -1 -31 -1 -32 -27 -15 -10 -44 -1 -13 -10 -24 -1 -21 -1 -32 -8 -9 -1 -16 -1 -7 -6 -12 -1 -5 -1 -4 -1 -1 0 0 0 1 1 4 1 5 1 12 6 7 1 16 1 9 8 32 1 21 1 24 10 13 1 44 10 15 27 32 1 31 1 80 14 19 12 60 1 21 16 68 1 41 1 48 39 25 1 112 14 45 20 56 1 81 16 92 22 31 1 92 1 33 51 192 18 61 1 72 26 59 1 156 1 39 55 80 18 71 1 176 108 43 1 124 22 45 32 140 1 123 20 96 34 49 24 272 1 77 75 140 D(10**1 ) / 7 == 1 D(10**2 ) / 7 == 20 D(10**3 ) / 7 == 300 D(10**4 ) / 7 == 4000 D(10**5 ) / 7 == 50000 D(10**6 ) / 7 == 600000 D(10**7 ) / 7 == 7000000 D(10**8 ) / 7 == 80000000 D(10**9 ) / 7 == 900000000 D(10**10) / 7 == 10000000000 D(10**11) / 7 == 110000000000 D(10**12) / 7 == 1200000000000 D(10**13) / 7 == 13000000000000 D(10**14) / 7 == 140000000000000 D(10**15) / 7 == 1500000000000000 D(10**16) / 7 == 16000000000000000 D(10**17) / 7 == 170000000000000000 D(10**18) / 7 == 1800000000000000000 D(10**19) / 7 == 19000000000000000000 D(10**20) / 7 == 200000000000000000000
Wren
As integer arithmetic in Wren is inaccurate above 2^53 we need to use BigInt here.
import "./big" for BigInt
import "./fmt" for Fmt
var D = Fn.new { |n|
if (n < 0) return -D.call(-n)
if (n < 2) return BigInt.zero
var f = BigInt.primeFactors(n)
var c = f.count
if (c == 1) return BigInt.one
if (c == 2) return f[0] + f[1]
var d = n / f[0]
return D.call(d) * f[0] + d
}
var ad = List.filled(200, 0)
for (n in -99..100) ad[n+99] = D.call(BigInt.new(n))
Fmt.tprint("$4i", ad, 10)
System.print()
for (m in 1..20) {
Fmt.print("D(10^$-2d) / 7 = $i", m, D.call(BigInt.ten.pow(m))/7)
}
- Output:
-75 -77 -1 -272 -24 -49 -34 -96 -20 -123 -1 -140 -32 -45 -22 -124 -1 -43 -108 -176 -1 -71 -18 -80 -55 -39 -1 -156 -1 -59 -26 -72 -1 -61 -18 -192 -51 -33 -1 -92 -1 -31 -22 -92 -16 -81 -1 -56 -20 -45 -14 -112 -1 -25 -39 -48 -1 -41 -1 -68 -16 -21 -1 -60 -12 -19 -14 -80 -1 -31 -1 -32 -27 -15 -10 -44 -1 -13 -10 -24 -1 -21 -1 -32 -8 -9 -1 -16 -1 -7 -6 -12 -1 -5 -1 -4 -1 -1 0 0 0 1 1 4 1 5 1 12 6 7 1 16 1 9 8 32 1 21 1 24 10 13 1 44 10 15 27 32 1 31 1 80 14 19 12 60 1 21 16 68 1 41 1 48 39 25 1 112 14 45 20 56 1 81 16 92 22 31 1 92 1 33 51 192 18 61 1 72 26 59 1 156 1 39 55 80 18 71 1 176 108 43 1 124 22 45 32 140 1 123 20 96 34 49 24 272 1 77 75 140 D(10^1 ) / 7 = 1 D(10^2 ) / 7 = 20 D(10^3 ) / 7 = 300 D(10^4 ) / 7 = 4000 D(10^5 ) / 7 = 50000 D(10^6 ) / 7 = 600000 D(10^7 ) / 7 = 7000000 D(10^8 ) / 7 = 80000000 D(10^9 ) / 7 = 900000000 D(10^10) / 7 = 10000000000 D(10^11) / 7 = 110000000000 D(10^12) / 7 = 1200000000000 D(10^13) / 7 = 13000000000000 D(10^14) / 7 = 140000000000000 D(10^15) / 7 = 1500000000000000 D(10^16) / 7 = 16000000000000000 D(10^17) / 7 = 170000000000000000 D(10^18) / 7 = 1800000000000000000 D(10^19) / 7 = 19000000000000000000 D(10^20) / 7 = 200000000000000000000