Évaluer la fonction Riemann Zeta à un nombre complexe

11

introduction

J'ai trouvé cette question qui était fermée parce qu'elle n'était pas claire, mais c'était une bonne idée. Je ferai de mon mieux pour en faire un défi clair.

La fonction Riemann Zeta est une fonction spéciale qui est définie comme la continuation analytique de

entrez la description de l'image ici

au plan complexe. Il existe de nombreuses formules équivalentes pour cela, ce qui le rend intéressant pour le golf de code.

Défi

Écrivez un programme qui prend 2 flottants en entrée (la partie réelle et imaginaire d'un nombre complexe) et évalue la fonction Riemann Zeta à ce point.

Règles

  • Entrée et sortie via console OU fonction entrée et valeur de retour
  • Les nombres complexes construits ne sont pas autorisés, utilisez des flotteurs (nombre, double, ...)
  • Pas de fonctions mathématiques sauf + - * / pow loget des fonctions trigonométriques réelles (si vous souhaitez intégrer, utilisez la fonction gamma, ... vous devez inclure cette définition de fonctions dans le code)
  • Entrée: 2 flotteurs
  • Sortie: 2 flotteurs
  • Votre code doit contenir une valeur qui donne une précision théoriquement arbitraire lorsqu'elle est rendue arbitraire grande / petite
  • Le comportement à l'entrée 1 n'est pas important (c'est le seul pôle de cette fonction)

Le code le plus court en octets gagne!

Exemple d'entrée et de sortie

Contribution:

2, 0

Production:

1,6449340668482266, 0

Contribution:

1, 1

Production:

0,5821580597520037, -0,9268485643308071

Contribution:

-dix

Production:

-0,08333333333333559, 0

Jens Renders
la source
1
Quelle est la précision de sortie requise? Je ne suis pas sûr de comprendre Votre code doit contenir une valeur qui donne une précision théoriquement arbitraire lorsqu'il est rendu arbitraire grand / petit . Voulez-vous dire comme une valeur maximale de boucle que lorsqu'elle est augmentée sans limite donne une précision accrue? Cette valeur peut-elle être codée en dur?
Luis Mendo
@DonMuesli Cela signifie que la précision dépend d'un paramètre, par exemple N, auquel vous pouvez donner la valeur que vous voulez, mais pour une précision donnée, vous pouvez rendre N petit ou assez grand pour atteindre cette précision. Le mot est théoriquement là parce que vous ne devez pas vous soucier de la précision limitée de la machine ou du langage.
Jens rend
Pour clarifier davantage N: est-il suffisant que pour toute borne epset entrée xexiste un Nqui calcule zeta(x)à l'intérieur eps; ou doit-il exister un Nqui ne dépend que de epset garantit que pour tout x(ou peut-être pour xplus qu'une fonction donnée epsdu pôle) il atteint la limite; ou peut Ndépendre x, mais les réponses devraient expliquer comment calculer Ndonné xet eps? (Ma théorie analytique des nombres n'est pas à la hauteur, mais je soupçonne que les options 2 et 3 vont dépasser toutes les affiches ordinaires sauf une ou deux).
Peter Taylor
@PeterTaylor N assez grand: pour tout xet pour tout, epsil doit exister un Ptel que pour toute N>Pla sortie soit plus proche que epsla valeur exacte. Est-ce clair? Dois-je le clarifier pour le cas avec N assez petit?
Jens Renders
Non, c'est assez clair.
Peter Taylor

Réponses:

8

Python - 385

Il s'agit d'une implémentation simple de l'équation 21 de http://mathworld.wolfram.com/RiemannZetaFunction.html Elle utilise la convention de Python pour les arguments facultatifs; si vous souhaitez spécifier une précision, vous pouvez passer un troisième argument à la fonction, sinon elle utilise 1e-24 par défaut.

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
la source
z(2,0)donne une valeur incorrecte, devrait être pi ^ 2/6.
GuillaumeDufay
4

Python 3 , 303 297 octets

Cette réponse est basée sur la réponse Python de RT avec plusieurs modifications:

  • Tout d'abord, Binomial(n, k)est défini comme celui p = p * (n-k) / (k+1)qui change Binomial(n,k)à Binomial(n,k+1)chaque passage de la boucle for.
  • Deuxièmement, (-1)**k * Binomial(n,k)est devenu p = p * (k-n) / (k+1)ce qui inverse le signe à chaque étape de la boucle for.
  • Troisièmement, la whileboucle a été modifiée pour vérifier immédiatement si a*a + b*b < E.
  • Quatrièmement, l'opérateur ne bitwise ~est utilisé dans plusieurs endroits où ils contribueraient à jouer au golf, en utilisant des identités telles que -n-1 == ~n, n+1 == -~net n-1 == ~-n.

Plusieurs autres petites modifications ont été apportées pour un meilleur golf, comme mettre la forboucle sur une ligne et l'appel printsur une ligne avec le code avant.

Suggestions de golf bienvenues. Essayez-le en ligne!

Edit: -6 octets à partir d'un certain nombre de petits changements.

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
la source
1

Axiom, 413 315 292 octets

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])

Cela implémenterait également l'équation 21 de http://mathworld.wolfram.com/RiemannZetaFunction.html Ce qui précède devrait être la fonction Axiom z (a, b) ici répétée 16 fois plus lentement que cette fonction ci-dessous Zeta (a, b) [ ce devrait être celui compilé] tout non golfé et commenté [1 seconde pour Zeta () contre 16 secondes pour z () pour une valeur de 20 chiffres après le point flottant]. Pour la question des chiffres, on choisirait la précision en appelant digits (); , par exemple les chiffres (10); z (1,1) devrait imprimer 10 chiffres après le point, mais les chiffres (50); z (1,1) devrait imprimer 50 chiffres après le point.

-- 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
la source