momente şi schiţe de informatică şi matematică
To attain knowledge, write. To attain wisdom, rewrite.

Asupra numărului de cifre, sub factoriale (IV)

R | factorial
2024 apr

[1] Asupra numărului de cifre, sub factoriale (III)

[2] Kevin Ryde - program C adăugat în martie 2024 la A084420

Reprezentarea în baza $10$ a numărului natural $h$ conține ncz(h)=floor(log10(h))+1=$\boldsymbol{\lfloor}\lg h\boldsymbol{\rfloor}+1$ cifre; amintim aceste proprietăți ale funcției floor$(\cdot)$=$\boldsymbol{\lfloor}\cdot\boldsymbol{\rfloor}$ („partea întreagă”):
$(i)\quad x-1 < \boldsymbol{\lfloor}x\boldsymbol{\rfloor} \le x$
$(ii)\quad \boldsymbol{\lfloor}x\boldsymbol{\rfloor} + \boldsymbol{\lfloor}y\boldsymbol{\rfloor} \le \boldsymbol{\lfloor}x+y\boldsymbol{\rfloor} \le \boldsymbol{\lfloor}x\boldsymbol{\rfloor} + \boldsymbol{\lfloor}y\boldsymbol{\rfloor}+1$

Pentru un număr natural $\boldsymbol{n}\ge 1$, fie $H$ unul dintre numerele naturale pentru care factorialul numărului $H+1$ are cu $\boldsymbol{n}$ cifre mai mult ca factorialul lui $H$: ncz$((H+1)!)$ = ncz$(H!)$ + $\boldsymbol{n}$; nu reușim să „caracterizăm” altfel, aceste numere – dar vom depista extremele.

Obs. Să remarcăm că în cazul a doi factori, chestiunea numărului de cifre ale produsului este banală: pentru $A$ și $B$ numere naturale cu $n_a$ și respectiv $n_b$ cifre, adică $10^{n_a-1}\le A < 10^{n_a}$ și $10^{n_b-1}\le B < 10^{n_b}$, avem $10^{n_a+n_b-2}\le A\times B < 10^{n_a+n_b}$, adică $A\times B$ are sau $n_a+n_b-1$, sau $n_a+n_b$ cifre (iar „cazurile extreme” sunt ușor de evidențiat).

Fiindcă $(H+1)!=H!\,(H+1)$, rezultă că $\boldsymbol{\lfloor}\lg(H!)+\lg(H+1)\boldsymbol{\rfloor}=\boldsymbol{\lfloor}\lg(H!)\boldsymbol{\rfloor}+n$ și ținând cont de $(ii)$ rezultă $\boldsymbol{\lfloor}\lg(H+1)\boldsymbol{\rfloor}\le n \le \boldsymbol{\lfloor}\lg(H+1)\boldsymbol{\rfloor}+1$; ținând cont de $(i)$ (cu $x=\lg(H+1)$), rezultă că $\lg(H+1)-1 < n \le \lg(H+1)+1$, adică $H+1<10^{n+1}$ și $H+1\ge 10^{n-1}$; totuși $H+1$ nu poate fi și $10^{n-1}$, fiindcă atunci creșterea de la $H!$ la $(H+1)!$ ar fi nu de $n$, ci numai de $n-1$ cifre: $\lg(10^{n-1}!)=\lg((10^{n-1}-1)!)+$ $\lg(10^{n-1})=\lg((10^{n-1}-1)!)+n\boldsymbol{-1}$.

Deci, dacă $(H+1)!$ are cu $\boldsymbol{n}$ cifre mai mult ca $H!$, atunci $\boldsymbol{10^{n-1}-1 < H < 10^{n+1}-1}$.

Intervalul $\mathcal{I}_n=\boldsymbol{(}10^{n-1}-1, 10^{n+1}-1\boldsymbol{)}$ conține $99\times 10^{n-1}-1$ numere naturale; cele aflate în intervalul închis $[10^{n-1}, 10^n-1]$ au câte $n$ cifre, iar următoarele până la $10^{n+1}-1$ au câte $n+1$ cifre. Deci pentru $h\in\mathcal{I}$, produsul $(h+1)!=h!\times(h+1)$ are fie cu $n-1$, fie cu $n$, fie cu $n+1$ cifre mai mult ca factorul $h!$.
Mai precis, cel puțin pentru câteva numere $h$ de la începutul intervalului $\mathcal{I}_n$, trecerea de la $h$ la $h+1$ (indicată mai jos prin "++") crește numărul cifrelor factorialului cu $n-1$, iar pentru câteva numere consecutive de la sfârșit, creșterea numărului de cifre ale factorialului devine egală cu $\boldsymbol{n+1}$:

$10^{n-1}-1$$++$$\cdots$$++$$H_1$$\boldsymbol{++}\,\cdots\quad\cdots\,\boldsymbol{++}$$H_2$$++$$\cdots$$++$$10^{n+1}-1$
$n-1$$\cdots$$n-1$$\boldsymbol{n}$$\quad n-1$ | $n$ | $n+1\quad$
(intercalat)
$\boldsymbol{n+1}$$n+1$$\cdots$$n+1$

Plecând de la capătul stâng $10^{n-1}-1$, după un anumit număr de incrementări succesive ajungem la primul (cel mai mic) număr $\boldsymbol{H_1}$ pentru care diferența numărului de cifre între $(H_1+1)!$ și $H_1!$ este egală cu $\boldsymbol{n}$ (și nu cu $n-1$, ca la pașii precedenți).
Mai precis, plecând de la $Z=10^{n-1}$, fie $k$ numărul de incrementări succesive (prin care factorialul crește cu câte $n-1$ cifre) după care ajungem la $H_1$:

    ncz((Z+1)!) = ncz(Z!) + n-1
    ncz((Z+2)!) = ncz((Z+1)!) + n-1
        ...
    ncz((Z+k)!) = ncz((Z+k-1)!) + n-1
   -----------------------------------  # ncz((Z+k)!) = ncz(Z!) + k(n-1)
ncz((Z+k+1)!) = ncz((Z+k)!) + n  # H1 = Z+k = 10^(n-1) + k 

Rezultă (în funcție de $k$) $H_1$ și – însumând primele $k$ egalități – $\mathrm{ncz}(H_1!)$:

$\quad(1)\quad$ $\boldsymbol{H_1=10^{n-1}+k}$  și  $\boldsymbol{\mathrm{ncz}(H_1!)=\mathrm{ncz}(10^{n-1}!)+k(n-1)}$

Analog, decrementând succesiv de la capătul din dreapta $10^{n+1}-1$, după un anumit număr $q$ de pași înapoi, ajungem la cel mai mare $H_2$ pentru care diferența numărului de cifre între $(H_2+1)!$ și $H_2!$ este egală cu $\boldsymbol{n}$ (și nu cu $n+1$, ca pe lanțul celor $q$ decrementări succesive). Pentru $H_2$ și pentru numărul cifrelor factorialului său avem (în funcție de $q$):

$\quad(2)\quad$ $\boldsymbol{H_2=10^{n+1}-1-q}$  și  $\boldsymbol{\mathrm{ncz}(H_2!)=\mathrm{ncz}((10^{n+1}-1)!)-q(n+1)+1}$

Nu avem idee cum s-ar identifica exact locurile $H_1$ și $H_2$ (altfel decât ca în [1], cercetând pas cu pas din capătul stâng și din cel drept); dar probabil că pot fi aproximate cumva…

Să observăm că pentru $h < H_1$, trecerea "++" de la $h$ la $h+1$ crește numai partea fracționară a diferenței dintre $\lg((h+1)!)$ și $\lg(h!)$ (partea întreagă fiind mereu $n-1$).
Deci pentru $h < H_1$ avem $\lg((h+1)!)-\lg(h!)=\lg(h+1) = (n-1)+\delta_h$, unde $0<\delta_h<1$ și $\delta_h$ crește odată cu $h$; pornind de la partea fracționară a factorialului de $10^{n-1}$, când $\sum\delta_h$ ajunge să depășească $1$ – atingem $h=H_1$ (fiindcă atunci partea întreagă urcă de la $n-1$ la $n$). Putem ilustra lucrurile cam așa, printr-un program $\boldsymbol{\mathcal{R}}$:

drafts <- function(n = 2) {
    s1 <- 10^(n-1); s2 <- 10^(n+1)-1
    H <- s1 : s2  # secvența de 99*10^(n-1)-1 numere (h+1)
    lgH <- log10(H)  # lg(H+1) ( = lg((H+1)!) - lg(H!) )
    d_ncz <- floor(lfactorial(H)/log(10)) - 
             floor(lfactorial(H-1)/log(10))  # ncz((H+1)!) - ncz(H!)
    deltaH <- lgH - d_ncz  # {lg(H+1) - ([lg((H+1)!)] - [lg(H!)])}
    for(s in c(s1, s1^2))  # delta_h pentru lg(10^(n-1)!) și lg((10^n)!)
        deltaH[s-s1+1] <- lfactorial(s)/log(10) - 
                              floor(lfactorial(s)/log(10))
    cbind(H, lgH, d_ncz, deltaH) |> as.data.frame()
}

D2 <- drafts()  # n = 2
print(rbind(D2[1:6, ], D2[91:95, ]))
         H      lgH d_ncz       deltaH
    1   10 1.000000     1  0.559763033  # lfactorial(10)/log(10) = 6.559763
    2   11 1.041393     1  0.041392685
    3   12 1.079181     1  0.079181246
    4   13 1.113943     1  0.113943352
    5   14 1.146128     1  0.146128036  # sum(D2$deltaH[1:5]) = 0.9404084
    6   15 1.176091     2 -0.823908741  # delta_h devine < 0
    91 100 2.000000     2  0.970003655  # lfactorial(100)/log(10) = 157.97
    92 101 2.004321     2  0.004321374
    93 102 2.008600     2  0.008600172
    94 103 2.012837     2  0.012837225  # sum(D2$deltaH[91:94]) = 0.9957624
    95 104 2.017033     3 -0.982966661  # delta_h devine < 0

Am selectat din D2 numai liniile care ne permit să deducem că 14 (și respectiv 103) este primul (respectiv, ultimul) $h$ pentru care $(h+1)!$ are cu două cifre (și respectiv, cu trei cifre) mai mult ca $h!$.

Avem $\delta=\{\lg(10!)\}+\sum_{h=11}^{h=14}\{\lg(h)\}\approx 0.94$ (unde "$\{.\}$" desemnează „partea fracționară”), iar apoi, $\delta + \{\lg(15)\}\approx \delta+0.176$ depășește $\boldsymbol{1}$ – încât, ajungând la $h=15$, partea întreagă devine 2.
La fel, $\{\lg(100!)\}+\sum_{h=101}^{h=103}\lg(h) + \{\lg(104)\}$ depășește $\boldsymbol{1}$ și la $h=104$ partea întreagă va urca de la 2 la 3 (și rămâne 3 până la capătul din dreapta al intervalului $\mathcal{I}_2$).

Sintetizând, am avea această caracterizare „teoretică”, pentru indexul $k$ al lui $H_1$ din (1)): $$k = \min\limits_{h}\left(\sum_{i=10^{n-1}}^{h}\boldsymbol{\{}\lg(i+1)\boldsymbol{\}} \boldsymbol{>} 1 - \boldsymbol{\{}\lg(10^{n-1}!)\boldsymbol{\}}\right) \tag{3}$$ ($H_1$ este cel mai mic $h$ pentru care suma părților fracționare ale logaritmilor de $i+1$ pentru $10^{n-1}\le i\le h$, depășește complementara părții fracționare a logaritmului factorialului de $10^{n-1}$)

Pentru $H_2$ (din (2)) putem proceda asemănător funcției drafts():

n <- 1  # =2, =3, etc.
q <- 10^(n+1) - 1
L <- lfactorial(q)/log(10) - floor(lfactorial(q)/log(10))
print(L)  # mantisa logaritmului factorialului de 10^(n+1)-1 
s <- 0
while(TRUE) {
    q <- q - 1
    s <- s + (n+1 - log10(q))  # însumează mantisele
    cat(q, " ", s, "\n")
    if(s > 1 - L) break  # până se complementează L (față de 1)
}
print(q-1); print(s)
[1] 0.9700037  # {lg((10^(n+1)-1)!)}
    98   0.008773924 
    97   0.02200219 
    96   0.03973096 
[1] 95  # H2 (la q=3), pentru n=1
[1] 0.03973096

Pentru $n>1$, coborârea pas cu pas până la $H_2$ este mult mai lungă:

[1] 0.4542745208  # {lg((10^(n+1)-1)!)}, pentru n=3
    9998   0.00008686758343 
    9997   0.0002171754752 
    9996   0.0003909280207 
        ...
    9843   0.5414624379 
    9842   0.5483790772 
[1] 9841  # H2 (la q=157), pentru n=3
[1] 0.5483790772

„Sintetizând” iarăși, am avea această formulă: $$H_2 = \max\limits_{q}\left(q(n+1) - \sum_{i=1}^{q}\boldsymbol{\{}\lg((10^{n+1}-1)-i)\boldsymbol{\}} \boldsymbol{>} 1 - \boldsymbol{\{}\lg((10^{n+1}-1)!)\boldsymbol{\}}\right) \tag{4}$$

$h$ din (3) și $q$ din (4) sunt mici față de capetele intervalului $\mathcal{I}_n$ – încât le vom putea aproxima, ținând seama că pentru $\boldsymbol{|x|<1}$ avem $\boldsymbol{\ln(1+x)}= x-x^2/2+x^3/3-\cdots \boldsymbol{\approx x}$ (v. [2], de unde am preluat aici ideea de estimare – în locul celei de căutare, din [1] – a valorilor $\boldsymbol{H_1}$ (după o discuție în cursul lunii februarie, pe care n-am reușit să o susținem până la capăt – dar iată că încercăm să o valorificăm aici – privitoare la secvența valorilor $\boldsymbol{H_2}$).

Notând $S_1=10^{n-1}$, suma părților fracționare din (3) se rescrie succesiv astfel: $\sum_{h=1}^{k}\boldsymbol{\{}\lg(S_1+h)\boldsymbol{\}}=\sum_{h=1}^{k}\boldsymbol{(}\lg(S_1+h)-\lfloor\,\lg(S_1+h)\,\rfloor\boldsymbol{)} = \sum_{h=1}^{k}\lg(S_1+h)-k(n-1)=$ $\sum_{h=1}^{k}\lg(S_1(1+h/S_1))-k(n-1)=\sum_{h=1}^{k}\lg S_1 + \sum_{h=1}^{k}\lg(1+h/S_1) - k(n-1)=$ $\sum_{h=1}^{k}\lg(1+h/S_1)=\sum_{h=1}^{k}(\ln(1+h/S_1)/\ln(10))$.
Aproximând $\ln(1+h/S_1)\approx h/S_1$, (3) revine acum la $\dfrac{1}{S_1}\dfrac{k(k+1)}{2} \approx (1-\delta_n)\ln(10)$, unde $\delta_n=\{\,\lg(10^{n-1}!)\,\}$; rădăcina pozitivă a ecuației $k^2+k-2S_1(1-\delta_n)\ln(10)\approx 0$ ne va da o aproximare destul de bună a indexului față de $S_1$ la care se află valoarea $H_1$:

$$H_1 \approx S_1 + \left\lfloor\,\dfrac{-1+\sqrt{1+8S_1(1-\delta_n)\ln(10)}}{2}\,\right\rceil\quad \text{unde}\, \delta_n=\{\,\lg(10^{n-1}!)\,\}\tag{5}$$

Probăm estimarea (5) prin următorul program:

estimate_H1 <- function(n) {
    S1 <- 10^(n-1)
    L <- lfactorial(S1) / log(10); delta <- L - floor(L)
    D <- 1 + 8*S1*(1 - delta)*log(10)
    S1 + round((-1 + sqrt(D))/2)
}
H1 <- sapply(1:13, estimate_H1)
print(H1)
 [1]            3           14          103         1042        10158
 [6]       100502      1000617     10006508    100019088   1000004378
[11]  10000170794 100000443062 1000001982570

Confruntând cu rezultatele cunoscute (v. [2]), constatăm că estimările obținute aplicând mot à mot formula (5) sunt chiar exacte pentru $n\le 7$; pentru $n=8$, $10$ și $11$ diferă cu o unitate (la ultima cifră, boldată mai sus); pentru $n=12$ diferența este de $+98$, iar pentru $n=13$, diferența este de $+914$. Este de așteptat o îmbunătățire sensibilă a estimărilor respective și pentru $n > 10$, dacă în loc de modelul numeric standard asumat de $\boldsymbol{\mathcal{R}}$, vom implica un model numeric prin care să mărim precizia calculelor intermediare (asigurând operanzilor un număr suficient de mare de cifre exacte); în orice caz, având prin (5) o estimare (totuși, acceptabilă) pentru $H_1$ – probabil că vom putea găsi valoarea exactă, mai sus sau mai jos de cea estimată, folosind „metoda înjumătățirii” (v. [2], unde în ecuația de gradul doi din care am dedus (5) s-a neglijat termenul de gradul întâi $k$, obținând în loc de (5) o formulă mai simplă, dar cu estimări mai largi pentru $H_1$ – încât s-a impus în mod firesc, metoda înjumătățirii).

Cam la fel procedăm acum pentru $H_2$, interpretând (4).
Notăm $S_2=10^{n+1}-1$, $L=\lg(S_2!)$ și $\delta_n=L-\lfloor\,L\,\rfloor$. Avem succesiv:

$\sum_{q=1}^{k}((n+1)-\lg(S_2-q)) = k(n+1)-\sum_{q=1}^{k}\lg(S_2(1-q/S_2)) =$ $k(n+1)-k\lg(S_2)-\sum_{q=1}^{k}\ln(1-q/S_2)/\ln(10) \boldsymbol{\approx} k(n+1-\lg(S_2))-\frac{1}{\ln(10)}\sum_{q=1}^{k}(-q/S_2) =$ $k(n+1-\lg(S_2))+\dfrac{1}{S_2\ln(10)}\dfrac{k(k+1)}{2} \boldsymbol{\approx} 1-\delta_n$.

Notând $a=\dfrac{1}{2S_2\ln(10)}$ și $b=n+1-\lg(S_2)+a$, rezultă „ecuația” $ak^2+bk-(1-\delta_n)\boldsymbol{\approx} 0$, a cărei rădăcină pozitivă, rotunjită: $\boldsymbol{k=\left\lfloor\,\dfrac{-b+\sqrt{b^2+4a(1-\delta_n)}}{2}\,\right\rceil}$ este indexul față de $S_2$ al unei valori din preajma lui $H_2$ (sau poate, chiar indexul lui $H_2$, măsurând de la $S_2-1$).
Prin următorul program verificăm formula obținută:

estimate_H2 <- function(n) {
    S2 <- 10^(n+1) - 1
    L <- lfactorial(S2) / log(10); d <- L - floor(L)
    a <- 1 / (2*S2*log(10))
    b <- n+1 - log10(S2) + a
    D <- b^2 + 4*a*(1-d)  # discriminantul ecuației ak^2 + bk - (1-d) = 0
    k <- (-b + sqrt(D)) / (2*a)
    S2 - round(k) - 1
}
H2 <- sapply(1:10, estimate_H2)
print(H2)
 [1]          96         957        9841       99497      999382     9993491
 [7]    99980911   999995621  9999829205 99999557000

Am boldat iarăși, acele ultime cifre la care avem diferențe față de valorile atestate în [1]; pentru $n=1$, valoarea exactă este 95 ($-1$ față de estimare), iar pentru $n=10$ diferența ar fi de $-29$ de poziții (în rest, ne-au rezultat chiar valorile $H_2$).

Pentru finalizare, găsind valorile exacte $H_1$ și $H_2$ (pentru $n$ cât se poate, de mare), rămâne de montat metoda înjumătățirii plecând de la formulele de estimare găsite și probate mai sus și desigur, de adoptat un model de calcul „în precizie arbitrară” („rămâne de …”; în realitate, ambele chestiuni rămase sunt dintre cele „delicate”).

vezi Cărţile mele (de programare)

docerpro | Prev | Next