Avaliar a função Riemann Zeta em um número complexo

11

Introdução

Encontrei essa pergunta que foi encerrada porque não estava clara, mas era uma boa ideia. Farei o meu melhor para tornar isso um desafio claro.

A função Riemann Zeta é uma função especial definida como a continuação analítica do

insira a descrição da imagem aqui

para o plano complexo. Existem muitas fórmulas equivalentes, o que o torna interessante para o código de golfe.

Desafio

Escreva um programa que use 2 carros alegóricos como entrada (a parte real e imaginária de um número complexo) e avalie a função Riemann Zeta nesse ponto.

Regras

  • Entrada e saída via console OU função de entrada e valor de retorno
  • Números complexos incorporados não são permitidos, use flutuadores (número, duplo, ...)
  • Nenhuma função matemática, exceto funções + - * / pow logtrigonométricas com valor real (se você deseja integrar, use a função gama, ... você deve incluir esta definição de funções no código)
  • Entrada: 2 carros alegóricos
  • Saída: 2 carros alegóricos
  • Seu código deve conter um valor que ofereça precisão teoricamente arbitrária quando convertido em grande / pequeno arbitrário
  • O comportamento na entrada 1 não é importante (este é o único pólo desta função)

O menor código em bytes vence!

Exemplo de entrada e saída

Entrada:

2, 0

Resultado:

1.6449340668482266, 0

Entrada:

1, 1

Resultado:

0,5821580597520037, -0,9268485643308071

Entrada:

-1, 0

Resultado:

-0.08333333333333559, 0

Jens Renders
fonte
1
Qual é a precisão de saída necessária? Não sei se entendi Seu código deve conter um valor que fornece precisão teoricamente arbitrária quando convertido em grande / pequeno arbitrário . Você quer dizer um valor máximo de loop do que quando aumentado sem limite fornece maior precisão? Esse valor pode ser codificado?
Luis Mendo
@DonMuesli Isso significa que a precisão depende de um parâmetro, por exemplo, N, que você pode atribuir qualquer valor que desejar, mas, para qualquer precisão, você pode tornar N pequeno ou grande o suficiente para atingir essa precisão. A palavra teoricamente está lá porque você não deve se preocupar com a precisão limitada da máquina ou do idioma.
Jens Renders
Para esclarecer melhor N: é suficiente que, para qualquer limite epse entrada x, exista um Nque calcule zeta(x)para dentro eps; ou deve existir um Nque dependa apenas epse garanta que, para qualquer x(ou talvez para xmais do que uma determinada função do epspólo) atinja o limite; ou pode Ndepender x, mas as respostas devem explicar como calcular Ndados xe eps? (Minha teoria dos números analíticos não é muito boa, mas suspeito que as opções 2 e 3 vão além de todos, exceto um ou dois pôsteres regulares).
Peter Taylor
@ PeterTaylor N grande o suficiente: para qualquer um xe para qualquer eps, deve existir um Ptal que, para toda N>Pa saída, seja mais próximo do epsque o valor exato. Isso está claro? Preciso esclarecê-lo para o caso com N pequeno o suficiente?
Jens Renders
Não, isso é claro o suficiente.
Peter Taylor

Respostas:

8

Python - 385

Esta é uma implementação direta da Equação 21 em http://mathworld.wolfram.com/RiemannZetaFunction.html Isso usa a convenção do Python para argumentos opcionais; se você deseja especificar uma precisão, pode passar um terceiro argumento para a função, caso contrário, ele usa 1e-24 por padrão.

import numpy as N
def z(r,i,E=1e-24):
 R=0;I=0;n=0;
 while(True):
  a=0;b=0;m=2**(-n-1)
  for k in range(0,n+1):
   M=(-1)**k*N.product([x/(x-(n-k))for x in range(n-k+1,n+1)]);A=(k+1)**-r;t=-i*N.log(k+1);a+=M*A*N.cos(t);b+=M*A*N.sin(t)
  a*=m;b*=m;R+=a;I+=b;n+=1
  if a*a+b*b<E:break
 A=2**(1-r);t=-i*N.log(2);a=1-A*N.cos(t);b=-A*N.sin(t);d=a*a+b*b;a=a/d;b=-b/d
 print(R*a-I*b,R*b+I*a)
RT
fonte
z(2,0)fornece um valor incorreto, deve ser pi ^ 2/6.
GuillaumeDufay
4

Python 3 , 303 297 bytes

Esta resposta é baseada na resposta Python da RT com várias modificações:

  • Primeiro, Binomial(n, k)é definido como o p = p * (n-k) / (k+1)que muda Binomial(n,k)para Binomial(n,k+1)cada passagem do loop for.
  • Segundo, (-1)**k * Binomial(n,k)tornou-se p = p * (k-n) / (k+1)que vira o sinal a cada passo do loop for.
  • Terceiro, o whileloop foi alterado para verificar imediatamente se a*a + b*b < E.
  • Em quarto lugar, o bit a bit operador não ~é usado em vários lugares onde eles iriam ajudar no golfe, usando identidades tais como -n-1 == ~n, n+1 == -~n, e n-1 == ~-n.

Várias outras pequenas modificações foram feitas para melhorar o golfe, como colocar o forloop em uma linha e a chamada printem uma linha com o código anterior.

Sugestões de golfe são bem-vindas. Experimente online!

Editar: -6 bytes de várias pequenas alterações.

import math as N
def z(r,i,E=1e-40):
 R=I=n=0;a=b=1
 while a*a+b*b>E:
  a=b=0;p=1;m=2**~n
  for k in range(1,n+2):M=p/k**r;p*=(k-1-n)/k;t=-i*N.log(k);a+=M*N.cos(t);b+=M*N.sin(t)
  a*=m;b*=m;R+=a;I+=b;n+=1
 A=2**-~-r;t=-i*N.log(2);x=1-A*N.cos(t);y=A*N.sin(t);d=x*x+y*y;return(R*x-I*y)/d,(R*y+I*x)/d
Sherlock9
fonte
1

Axiom, 413 315 292 bytes

p(n,a,b)==(x:=log(n^b);y:=n^a;[y*cos(x),y*sin(x)]);z(a,b)==(r:=[0.,0.];e:=10^-digits();t:=p(2,1-a,-b);y:=(1-t.1)^2+t.2^2;y=0=>[];m:=(1-t.1)/y;q:=t.2/y;n:=0;repeat(w:=2^(-n-1);abs(w)<e=>break;r:=r+w*reduce(+,[(-1)^k*binomial(n,k)*p(k+1,-a,-b) for k in 0..n]);n:=n+1);[r.1*m-q*r.2,m*r.2+r.1*q])

Isso também implementaria a equação 21 de http://mathworld.wolfram.com/RiemannZetaFunction.html A descrição acima deve ser a função Axiom repetida z (a, b) aqui 16x mais lenta que esta abaixo da função Zeta (a, b) [ esse deve ser o que foi compilado] todos sem golfe e comentado [1 segundo para Zeta () contra 16 segundos para z () para um valor de 20 dígitos após o ponto de flutuação]. Para a pergunta do dígito, um escolheria a precisão chamando dígitos (); função, por exemplo, dígitos (10); z (1,1) deve imprimir 10 dígitos após o ponto, mas dígitos (50); z (1,1) devem imprimir 50 dígitos após o ponto.

-- elevImm(n,a,b)=n^(a+i*b)=r+i*v=[r,v]
elevImm(n:INT,a:Float,b:Float):Vector Float==(x:=log(n^b);y:=n^a;[y*cos(x),y*sin(x)]::Vector Float);

--                      +oo               n
--                      ---              ---
--             1        \       1        \            n 
--zeta(s)= ---------- * /     ------  *  /    (-1)^k(   )(k+1)^(-s)
--          1-2^(1-s)   ---n  2^(n+1)    ---k         k  
--                       0                0


Zeta(a:Float,b:Float):List Float==
  r:Vector Float:=[0.,0.]; e:=10^-digits()

  -- 1/(1-2^(1-s))=1/(1-x-i*y)=(1-x+iy)/((1-x)^2+y^2)=(1-x)/((1-x)^2+y^2)+i*y/((1-x)^2+y^2)    

  t:=elevImm(2,1-a,-b);
  y:=(1-t.1)^2+t.2^2;
  y=0=>[] 
  m:=(1-t.1)/y; 
  q:=t.2/y
  n:=0
  repeat
     w:=2^(-n-1)
     abs(w)<e=>break  --- this always terminate because n increase
     r:=r+w*reduce(+,[(-1)^k*binomial(n,k)*elevImm(k+1,-a,-b) for k in 0..n])
     n:=n+1
  -- (m+iq)(r1+ir2)=(m*r1-q*r2)+i(m*r2+q*r1)
  [r.1*m-q*r.2,m*r.2+r.1*q]

this is one test for the z(a,b) function above:

(10) -> z(2,0)
   (10)  [1.6449340668 482264365,0.0]
                                              Type: List Expression Float
(11) -> z(1,1)
   (11)  [0.5821580597 520036482,- 0.9268485643 3080707654]
                                              Type: List Expression Float
(12) -> z(-1,0)
   (12)  [- 0.0833333333 3333333333 3,0.0]
                                              Type: List Expression Float
(13) -> z(1,0)
   (13)  []
RosLuP
fonte