Execuções de dígitos no Pi

13

Seu objetivo é produzir a sequência estritamente crescente de dígitos consecutivos idênticos de pi (π). Cada termo na sequência deve ter um dígito a mais que o anterior. Então 3(0º dígito de pi) é a primeira vez que ocorre uma sequência de dígitos (comprimento 1). O próximo a ocorrer é 33(dígitos 24 e 25 de pi). Obviamente, essa sequência requer que os dígitos de pi estejam na base 10 .

Os conhecidos até agora e os seis primeiros ocorrem nos primeiros 800 dígitos:

3
33
111
9999
99999
999999
3333333
44444444
777777777
6666666666
... (not in first 2 billion digits)

Observe que todos os nove consecutivos ocorrem juntos, na mesma execução, portanto, se a próxima execução maior que você encontrou fosse 1000 0s consecutivos , isso preencheria vários termos da sequência.

Não encontrei mais termos com o meu programa. Sei que não há mais termos nos primeiros 50000 dígitos ou mais. Meu programa estava demorando muito com 500000 dígitos, então desisti.

Implementação de referência

Você pode:

  • Saída da sequência para sempre
  • Pegue um número inteiro ne encontre os primeiros nnúmeros na sequência
  • Pegue um número inteiro ne encontre os números na sequência contida nos primeiros ndígitos de pi.

Certifique-se de especificar qual é o seu código. O número npode ser zero ou um indexado.

Inspirado por esta pergunta do mathoverflow .

mbomb007
fonte
1
Relacionado - que séries de 9s causaram dores de cabeça por muitas respostas: P
Mego
Você tem permissão para iniciar a saída com a sequência vazia?
LegionMammal978
2
Além disso, o próximo termo da sequência parece ser 3333333 nos dígitos 10 ^ -710100 a 10 ^ -710106. O valor para n = 8 não aparece nos primeiros 5 000 000 dígitos.
LegionMammal978
4
Mais dois termos: 44444444 nos dígitos 10 ^ -22931745 a 10 ^ -22931752 e 777777777 nos dígitos 10 ^ -24658601 a 10 ^ -24658609. O valor para n = 10 não aparece nos primeiros 100 000 000 dígitos.
LegionMammal978
1
Mais um termo: 6666666666 em 10 ^ -386980412. O décimo primeiro termo não aparece nos primeiros 2 000 000 000 dígitos.
Primo

Respostas:

5

Mathematica, 85 bytes

FromDigits/@DeleteDuplicatesBy[Join@@Subsets/@Split@RealDigits[Pi,10,#][[1]],Length]&

Função anônima. Pega n como entrada e retorna os elementos da sequência nos primeiros n dígitos de π. A saída está na forma de {0, 3, 33, 111, ...}.

LegionMammal978
fonte
4

Python 2, 110 bytes

n=input()
x=p=7*n|1
while~-p:x=p/2*x/p+2*10**n;p-=2
l=m=0
for c in`x`:
 l=l*(p==c)+1;p=c
 if l>m:m=l;print p*l

O número máximo de dígitos a serem verificados é obtido a partir do stdin. 10.000 dígitos terminam em cerca de 2s com o PyPy 5.3.

Uso da amostra

$ echo 10000 | pypy pi-runs.py
3
33
111
9999
99999
999999

Algo util

from sys import argv
from gmpy2 import mpz

def pibs(a, b):
  if a == b:
    if a == 0:
      return (1, 1, 1123)
    p = a*(a*(32*a-48)+22)-3
    q = a*a*a*24893568
    t = 21460*a+1123
    return (p, -q, p*t)
  m = (a+b) >> 1
  p1, q1, t1 = pibs(a, m)
  p2, q2, t2 = pibs(m+1, b)
  return (p1*p2, q1*q2, q2*t1 + p1*t2)

if __name__ == '__main__':
  from sys import argv
  digits = int(argv[1])

  pi_terms = mpz(digits*0.16975227728583067)
  p, q, t = pibs(0, pi_terms)

  z = mpz(10)**digits
  pi = 3528*q*z/t

  l=m=0
  x=0
  for c in str(pi):
   l=l*(p==c)+1;p=c
   if l>m:m=l;print x,p*l
   x+=1

Eu mudei de Chudnovsky para Ramanujan 39 para isso. Chudnovsky ficou sem memória no meu sistema logo após 100 milhões de dígitos, mas Ramanujan chegou a 400 milhões, em apenas 38 minutos. Acho que esse é outro caso em que a taxa de crescimento mais lenta dos termos vence no final, pelo menos em um sistema com recursos limitados.

Uso da amostra

$ python pi-ramanujan39-runs.py 400000000
0 3
25 33
155 111
765 9999
766 99999
767 999999
710106 3333333
22931752 44444444
24658609 777777777
386980421 6666666666

Geradores não ligados mais rápidos

A implementação de referência dada na descrição do problema é interessante. Ele usa um gerador ilimitado, extraído diretamente do algoritmo de torneira não encadernada para os dígitos do Pi . Segundo o autor, as implementações fornecidas são "deliberadamente obscuras", por isso decidi fazer novas implementações dos três algoritmos listados pelo autor, sem ofuscação deliberada. Eu também adicionei um quarto, baseado no Ramanujan # 39 .

try:
  from gmpy2 import mpz
except:
  mpz = long

def g1_ref():
  # Leibniz/Euler, reference
  q, r, t = mpz(1), mpz(0), mpz(1)
  i, j = 1, 3
  while True:
    n = (q+r)/t
    if n*t > 4*q+r-t:
      yield n
      q, r = 10*q, 10*(r-n*t)
    q, r, t = q*i, (2*q+r)*j, t*j
    i += 1; j += 2

def g1_md():
  # Leibniz/Euler, multi-digit
  q, r, t = mpz(1), mpz(0), mpz(1)
  i, j = 1, 3
  z = mpz(10)**10
  while True:
    n = (q+r)/t
    if n*t > 4*q+r-t:
      for d in digits(n, i>34 and 10 or 1): yield d
      q, r = z*q, z*(r-n*t)
    u, v, x = 1, 0, 1
    for k in range(33):
      u, v, x = u*i, (2*u+v)*j, x*j
      i += 1; j += 2
    q, r, t = q*u, q*v+r*x, t*x

def g2_md():
  # Lambert, multi-digit
  q, r, s, t = mpz(0), mpz(4), mpz(1), mpz(0)
  i, j, k = 1, 1, 1
  z = mpz(10)**49
  while True:
    n = (q+r)/(s+t)
    if n == q/s:
      for d in digits(n, i>65 and 49 or 1): yield d
      q, r = z*(q-n*s), z*(r-n*t)
    u, v, w, x = 1, 0, 0, 1
    for l in range(64):
      u, v, w, x = u*j+v, u*k, w*j+x, w*k
      i += 1; j += 2; k += j
    q, r, s, t = q*u+r*w, q*v+r*x, s*u+t*w, s*v+t*x

def g3_ref():
  # Gosper, reference
  q, r, t = mpz(1), mpz(180), mpz(60)
  i = 2
  while True:
    u, y = i*(i*27+27)+6, (q+r)/t
    yield y
    q, r, t, i = 10*q*i*(2*i-1), 10*u*(q*(5*i-2)+r-y*t), t*u, i+1

def g3_md():
  # Gosper, multi-digit
  q, r, t = mpz(1), mpz(0), mpz(1)
  i, j = 1, 60
  z = mpz(10)**50
  while True:
    n = (q+r)/t
    if n*t > 6*i*q+r-t:
      for d in digits(n, i>38 and 50 or 1): yield d
      q, r = z*q, z*(r-n*t)
    u, v, x = 1, 0, 1
    for k in range(37):
      u, v, x = u*i*(2*i-1), j*(u*(5*i-2)+v), x*j
      i += 1; j += 54*i
    q, r, t = q*u, q*v+r*x, t*x

def g4_md():
  # Ramanujan 39, multi-digit
  q, r, s ,t = mpz(0), mpz(3528), mpz(1), mpz(0)
  i = 1
  z = mpz(10)**3511
  while True:
    n = (q+r)/(s+t)
    if n == (22583*i*q+r)/(22583*i*s+t):
      for d in digits(n, i>597 and 3511 or 1): yield d
      q, r = z*(q-n*s), z*(r-n*t)
    u, v, x = mpz(1), mpz(0), mpz(1)
    for k in range(596):
      c, d, f = i*(i*(i*32-48)+22)-3, 21460*i-20337, -i*i*i*24893568
      u, v, x = u*c, (u*d+v)*f, x*f
      i += 1
    q, r, s, t = q*u, q*v+r*x, s*u, s*v+t*x

def digits(x, n):
  o = []
  for k in range(n):
    x, r = divmod(x, 10)
    o.append(r)
  return reversed(o)

Notas

Acima estão 6 implementações: as duas implementações de referência fornecidas pelo autor (denotadas _ref) e quatro que computam termos em lotes, gerando vários dígitos ao mesmo tempo ( _md). Todas as implementações foram confirmadas para 100.000 dígitos. Ao escolher tamanhos de lote, escolhi valores que lentamente perdem a precisão ao longo do tempo. Por exemplo, g1_mdgera 10 dígitos por lote, com 33 iterações. No entanto, isso produzirá apenas ~ 9,93 dígitos corretos. Quando a precisão acabar, a condição de verificação falhará, acionando um lote extra a ser executado. Isso parece ser mais eficiente do que a precisão desnecessariamente desnecessária e desnecessária ao longo do tempo.

  • g1 (Leibniz / Euler)
    Uma variável extra jé mantida, representando 2*i+1. O autor faz o mesmo na implementação de referência. Calcular nseparadamente é muito mais simples (e menos obscuro), porque usa os valores atuais de q, re t, em vez do seguinte.
  • g2 (Lambert)
    O cheque n == q/sé reconhecidamente bastante relaxado. Isso deve ler n == (q*(k+2*j+4)+r)/(s*(k+2*j+4)+t), onde jestá 2*i-1e kestá i*i. Em iterações mais altas, os termos re tse tornam cada vez menos significativos. Como é, é bom para os primeiros 100.000 dígitos, então provavelmente é bom para todos. O autor não fornece nenhuma implementação de referência.
  • g3 (Gosper)
    O autor conjectura que não é necessário verificar se isso nnão será alterado nas iterações subseqüentes e que serve apenas para retardar o algoritmo. Embora provavelmente seja verdade, o gerador está mantendo ~ 13% mais dígitos corretos do que o gerado atualmente, o que parece um pouco inútil. Adicionei o check-in novamente e aguarde até 50 dígitos estarem corretos, gerando todos de uma só vez, com um ganho notável no desempenho.
  • g4 (Ramanujan 39)
    Calculado como

    Infelizmente, snão é zerado, devido à composição inicial (3528 ÷), mas ainda é significativamente mais rápido que o g3. A convergência é de ~ 5,89 dígitos por termo, 3511 dígitos são gerados por vez. Se isso é um pouco demais, gerar 271 dígitos por 46 iterações também é uma escolha decente.

Horários

Tirada no meu sistema, apenas para fins de comparação. Os horários são listados em segundos. Se um tempo demorasse mais de 10 minutos, não executei mais testes.

            |  g1_ref |  g1_md  |  g2_md  |  g3_ref |  g3_md  |  g4_md 
------------+---------+---------+---------+---------+---------+--------
    10,000  |  1.645  |  0.229  |  0.093  |  0.312  |  0.062  |  0.062 
    20,000  |  6.859  |  0.937  |  0.234  |  1.140  |  0.250  |  0.109 
    50,000  |  55.62  |  5.546  |  1.437  |  9.703  |  1.468  |  0.234 
   100,000  |  247.9  |  24.42  |  5.812  |  39.32  |  5.765  |  0.593 
   200,000  |  2,158  |  158.7  |  25.73  |  174.5  |  33.62  |  2.156 
   500,000  |    -    |  1,270  |  215.5  |  3,173  |  874.8  |  13.51 
 1,000,000  |    -    |    -    |  1,019  |    -    |    -    |  58.02 

É interessante que, g2eventualmente g3, supere , apesar de uma taxa mais lenta de convergência. Suspeito que isso ocorra porque os operandos crescem a uma taxa significativamente mais lenta, vencendo a longo prazo. A implementação mais rápida g4_mdé aproximadamente 235x mais rápida que a g3_refimplementação em 500.000 dígitos. Dito isto, ainda há uma sobrecarga significativa para a transmissão de dígitos dessa maneira. O cálculo de todos os dígitos diretamente usando o Ramanujan 39 ( fonte python ) é aproximadamente 10 vezes mais rápido.

Por que não Chudnovsky?

O algoritmo de Chudnovsky requer uma raiz quadrada de precisão total, na qual sinceramente não sei como trabalhar - assumindo que possa ser. Ramanujan 39 é um tanto especial nesse sentido. No entanto, o método parece propício a fórmulas do tipo Machin, como as usadas pelo triturador de y, de modo que pode ser uma avenida que vale a pena explorar.

primo
fonte
TIL Ideone suporta Pypy. Então, o segundo programa é construído para velocidade?
mbomb007
@ mbomb007 "Então, o segundo programa foi desenvolvido para velocidade?" Isto é. Eu acho que o desafio teria sido tão interessante quanto um código mais rápido .
Primo
Mesmo. Eu considerei os dois. Não sei como as pessoas se sentem ao publicar novamente com uma tag diferente. Ele pode ser mais útil se for para ser adicionado ao OEIS (que não contém esta sequência)
mbomb007
3

Haskell, 231 bytes

import Data.List
g(q,r,t,k,n,l)|4*q+r-t<n*t=n:g(10*q,10*(r-n*t),t,k,div(10*(3*q+r))t-10*n,l)|0<1=g(q*k,(2*q+r)*l,t*l,k+1,div(q*(7*k+2)+r*l)(t*l),l+2)
p=nubBy(\x y->length x==length y).concatMap inits.group$g(1,0,1,1,3,3) 

Isso usa os algoritmos de espigão ilimitado para os dígitos do Pi de Jeremy Gibbons, 2004. O resultado é p. Tecnicamente, ele deve suportar infinitas seqüências de saída, mas isso pode demorar um pouco (e é limitado pela sua memória).

Zeta
fonte
3

Python 2, 298 bytes

Observe que o código para gerar pi é obtido da implementação do OP.

def p():
 q,r,t,j=1,180,60,2
 while 1:
  u,y=3*(3*j+1)*(3*j+2),(q*(27*j-12)+5*r)//(5*t)
  yield y
  q,r,t,j=10*q*j*(2*j-1),10*u*(q*(5*j-2)+r-y*t),t*u,j+1
p=p()
c=r=0
d=[0]
while 1:
 t=p.next()
 if t==d[len(d)-1]:d.append(t)
 else:d=[t]
 if len(d)>r:r=len(d);print"".join([`int(x)`for x in d])
 c+=1

Minha primeira tentativa de jogar golfe em Python. Produz a sequência para sempre.

acrólito
fonte
Poderia explicar como você calcula πaqui? Você, é claro, calcula pi, certo?
R. Kap 14/07
Não é possível testar agora, mas você não está calculando πpara sempre?
Yytsi 14/07/16
Não @TuukkaX não aparecer, de modo que tem um yieldque pára-lo, mas eu não sou muito bom em python
Downgoat
Downgoat está correto - ele usa uma função de gerador .
Mego
1
Eu escrevi todo o código, não olhei para sua implementação, exceto a pparte
acrolith 14/07/16
3

Python 3.5, 278 263 bytes:

import decimal,re;decimal.getcontext().prec=int(input());D=decimal.Decimal;a=p=1;b,t=1/D(2).sqrt(),1/D(4)
for i in[1]*50:z=(a+b)/2;b=(a*b).sqrt();t-=p*(a-z)**2;a=z;p*=2;pi=(z*2)**2/(4*t);i=0;C=lambda r:re.search(r'(\d)\1{%s}'%r,str(pi))
while C(i):print(C(i));i+=1

Isso recebe ncomo entrada os primeiros ndígitos de πe, em seguida, gera os membros da sequência nesses primeiros ndígitos. Agora, isso usa o módulo decimal embutido do Python para ir além das limitações de ponto flutuante do Python e, em seguida, define a precisão, ou epsilon, para o quanto as entradas do usuário. Então, para calcular π, isso passa por 50 iterações usando o eficiente algoritmo de Gausse-Legendre , já que o algoritmo aparentemente dobra o número de dígitos corretos a cada vez e, portanto, em 50 iterações, podemos obter 2^50ou 1,125,899,906,842,624corrigir dígitos. Finalmente, após os cálculos, ele usa uma expressão regular com formatação de string em um whileloop para encontrar e imprimirre objetos de correspondência (que espero que estejam bem) para todos os dígitos recorrentes e contínuos 1 dígito a mais do que na iteração anterior através do loop.

Consegui usar esse algoritmo para calcular com êxito e precisão πaté 10,000,000(dez milhões) dígitos, o que levou cerca de 4 horas e 12 minutos para ser concluído. O seguinte foi a saída final:

<_sre.SRE_Match object; span=(0, 1), match='3'>
<_sre.SRE_Match object; span=(25, 27), match='33'>
<_sre.SRE_Match object; span=(154, 157), match='111'>
<_sre.SRE_Match object; span=(763, 767), match='9999'>
<_sre.SRE_Match object; span=(763, 768), match='99999'>
<_sre.SRE_Match object; span=(763, 769), match='999999'>
<_sre.SRE_Match object; span=(710101, 710108), match='3333333'> 

Então, posso dizer com confiança que o 8º número da sequência nem ocorre nos primeiros 10 milhões de dígitos! πé um número aleatório ...

R. Kap
fonte