Implémenter la fonction Polygamma

Question

La fonction Polygamma d'ordre \$m\$, \$\psi^{(m)}(z)\$, est le \$(m + 1)\$ème dérivée du logarithme de la fonction gamma , qui est aussi le \$m\$ème dérivé de la fonction digamma . Votre tâche est de prendre un entier \$m\$et un nombre réel positif \$z\$et sortie \$\psi^{(m)}(z)\$

Définitions

Pour ceux qui ne connaissent pas les fonctions ci-dessus (Gamma, Digamma et Polygamma), voici quelques définitions différentes pour chacune:

\$\Gamma(z)\$

• La fonction gamma est une extension de la factorielle ( \$x! = 1\cdot2\cdot3\cdots(x-1)\cdot(x)\$) en nombres réels
• \$\Gamma(z) = \int_{0}^{\infty}x^{z-1}e^{-x}dx\$
• \$\Gamma(n) = (n - 1)! \:,\:\: n \in \mathbb{N}\$
• \$\Gamma(n+1) = n\Gamma(n) \:,\:\: n \in \mathbb{N}\$

\$\psi(z)\$

• La fonction digamma est le dérivé logarithmique de la fonction gamma
• \$\psi(z) = \frac{d}{dz}\ln(\Gamma(z))\$
• \$\psi(z) = \frac{\Gamma'(z)}{\Gamma(z)}\$
• \$\psi(z + 1) = \psi(z) + \frac{1}{z}\$

\$\psi^{(m)}(z)\$

• La fonction polygamma d'ordre \$m\$est le \$m\$ème dérivée de la fonction digamma
• \$\psi^{(m)}(z) = \frac{d^m}{dz^m}\psi(z)\$
• \$\psi^{(m)}(z) = \frac{d^{m+1}}{dz^{m+1}}\ln(\Gamma(z))\$
• \$\psi^{(m)}(z+1)= \psi^{(m)}(z) + (-1)^m\frac{m!}{z^{m+1}}\$

Tâche

Vous devez prendre deux entrées, un nombre naturel \$m\$et un nombre réel positif \$z\$, et sortie \$\psi^{(m)}(z)\$. Les entrées et sorties s'inscriront toujours dans les limites numériques de votre langue, mais votre algorithme doit théoriquement fonctionner pour toutes les entrées .

Comme la sortie est généralement un nombre réel, plutôt qu'un entier, la sortie doit être correcte à au moins 10 chiffres significatifs. Les zéros de fin peuvent être omis pour les valeurs exactes. Par exemple, si la sortie est un entier, les décimales de fin 0ne sont pas obligatoires, mais sont autorisées si vous le souhaitez.

C'est du code-golf donc le code le plus court en octets l'emporte.

Cas de test

Les résultats peuvent différer en raison d'inexactitudes en virgule flottante, la bibliothèque scipy de Python a été utilisée pour générer les valeurs. Les valeurs sont arrondies à 15d.p., sauf indication contraire.

 m,                  z -> ψ⁽ᵐ⁾(z)
17,                  2 -> 1357763223.715975761413574
 5,                 40 -> 0.0000002493894351
 9,           53.59375 -> 0.00000000001201026493
35,                  9 -> 469354.958166260155849
46,                  5 -> -7745723758939047727202304.000000000000000
 7, 1.2222222222222222 -> 1021.084176496877490
28,               6.25 -> -2567975.924144014250487
 2,               7.85 -> -0.018426049840992

Cette table a les valeurs de \$\psi^{(m)}(z)\$pour \$0 \le m \le 9\$et \$1 \le z \le 20\$:

+---+------------------------+---------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+
|   |           1            |          2          |         3          |         4          |         5          |         6          |         7          |         8          |         9          |         10         |         11         |         12         |         13         |         14         |         15         |         16         |         17         |         18         |         19         |         20         |
+---+------------------------+---------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+
| 0 |   -0.577215664901533   |  0.422784335098467  | 0.922784335098467  | 1.256117668431800  | 1.506117668431800  | 1.706117668431800  | 1.872784335098467  | 2.015641477955610  | 2.140641477955610  | 2.251752589066721  | 2.351752589066721  | 2.442661679975812  | 2.525995013309145  | 2.602918090232222  | 2.674346661660794  | 2.741013328327460  | 2.803513328327460  | 2.862336857739225  | 2.917892413294781  | 2.970523992242149  |
| 1 |   1.644934066848227    |  0.644934066848227  | 0.394934066848226  | 0.283822955737115  | 0.221322955737115  | 0.181322955737115  | 0.153545177959338  | 0.133137014694031  | 0.117512014694031  | 0.105166335681686  | 0.095166335681686  | 0.086901872871768  | 0.079957428427324  | 0.074040268664010  | 0.068938227847684  | 0.064493783403239  | 0.060587533403239  | 0.057127325790783  | 0.054040906037696  | 0.051270822935203  |
| 2 |   -2.404113806319188   |  -0.404113806319189 | -0.154113806319189 | -0.080039732245115 | -0.048789732245114 | -0.032789732245115 | -0.023530472985855 | -0.017699569195768 | -0.013793319195768 | -0.011049834970802 | -0.009049834970802 | -0.007547205368999 | -0.006389797961592 | -0.005479465690312 | -0.004750602716552 | -0.004158010123959 | -0.003669728873959 | -0.003262645625435 | -0.002919710097314 | -0.002628122402315 |
| 3 |   6.493939402266829    |  0.493939402266829  | 0.118939402266829  | 0.044865328192755  | 0.021427828192755  | 0.011827828192755  | 0.007198198563125  | 0.004699239795945  | 0.003234396045945  | 0.002319901304290  | 0.001719901304290  | 0.001310093231071  | 0.001020741379219  | 0.000810664701232  | 0.000654479778283  | 0.000535961259764  | 0.000444408525389  | 0.000372570305061  | 0.000315414383708  | 0.000269374221340  |
| 4 |  -24.886266123440890   |  -0.886266123440879 | -0.136266123440878 | -0.037500691342113 | -0.014063191342113 | -0.006383191342113 | -0.003296771589026 | -0.001868795150638 | -0.001136373275638 | -0.000729931168235 | -0.000489931168235 | -0.000340910050701 | -0.000244459433417 | -0.000179820455575 | -0.000135196191875 | -0.000103591253604 | -0.000080703070010 | -0.000063799959344 | -0.000051098643488 | -0.000041405977726 |
| 5 |  122.081167438133861   |  2.081167438133896  | 0.206167438133897  | 0.041558384635954  | 0.012261509635954  | 0.004581509635954  | 0.002009493175049  | 0.000989510004771  | 0.000531746332896  | 0.000305945162117  | 0.000185945162117  | 0.000118208290511  | 0.000078020533309  | 0.000053159387985  | 0.000037222150950  | 0.000026687171526  | 0.000019534614153  | 0.000014563111016  | 0.000011034967722  | 0.000008484266206  |
| 6 |  -726.011479714984489  |  -6.011479714984437 | -0.386479714984435 | -0.057261607988551 | -0.013316295488551 | -0.004100295488551 | -0.001528279027645 | -0.000654007738836 | -0.000310684984930 | -0.000160150871077 | -0.000088150871077 | -0.000051203486564 | -0.000031109607963 | -0.000019635233198 | -0.000012804988755 | -0.000008590996985 | -0.000005908787970 | -0.000004154139804 | -0.000002978092040 | -0.000002172607350 |
| 7 |  5060.549875237640663  |  20.549875237639476 | 0.862375237639470  | 0.094199654649073  | 0.017295357774073  | 0.004392957774073  | 0.001392271903016  | 0.000518000614207  | 0.000217593204539  | 0.000100511115987  | 0.000050111115987  | 0.000026599144024  | 0.000014877714841  | 0.000008699205352  | 0.000005284083130  | 0.000003317553637  | 0.000002144087193  | 0.000001421585007  | 0.000000964233099  | 0.000000667475582  |
| 8 | -40400.978398747647589 | -80.978398747634884 | -2.228398747634885 | -0.179930526327158 | -0.026121932577158 | -0.005478092577158 | -0.001477178082416 | -0.000478010895205 | -0.000177603485537 | -0.000073530517936 | -0.000033210517936 | -0.000016110901963 | -0.000008296615840 | -0.000004494456155 | -0.000002542957742 | -0.000001494142013 | -0.000000907408791 | -0.000000567407762 | -0.000000364140247 | -0.000000239189714 |
| 9 | 363240.911422382690944 | 360.911422382626938 | 6.536422382626807  | 0.391017718703625  | 0.044948382766125  | 0.007789470766125  | 0.001788099024012  | 0.000503455497598  | 0.000165497161722  | 0.000061424194120  | 0.000025136194120  | 0.000011145599233  | 0.000005284884641  | 0.000002652620244  | 0.000001398085550  | 0.000000768796112  | 0.000000438758675  | 0.000000258758130  | 0.000000157124373  | 0.000000097937278  |
+---+------------------------+---------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+

Answer 1

JavaScript (ES7),  68 66 61  59 octets

Attend (m)(z).

(m,n=m)=>g=z=>n?-n--*g(z):eval("for(k=5e6;k--;)n-=z++**~m")

Essayez-le en ligne!

Ceci est basé sur la représentation de série suivante ( de Wikipedia ):

$$\psi^{(m)}(z)=(-1)^{m+1}m!\sum_{k=0}^{\infty}\dfrac{1}{(z+k)^{m+1}}$$

Commenté

(m, n = m) =>                // outer function taking m and saving a copy in n
g = z =>                     // inner recursive function taking z
  n ?                        // if n is not equal to 0:
    -n--                     //   yield -n to invert the sign; decrement n afterwards
    * g(z)                   //   multiply by the result of a recursive call
  :                          // else:
    eval(                    //   evaluate as JS code:
      "for(k = 5e6; k--;)" + //     repeat 5 million times:
        "n -= z++ ** ~m"     //       subtract z ** -(m+1) from n; increment z
    )                        //   end of eval(), which returns the final value of n

Answer 2

Mathematica, 9 octets

PolyGamma

Bien sûr, Mathematica a une fonction intégrée pour cela

Answer 3

R , 8 octets

psigamma

Essayez-le en ligne!

Prend des entrées z,m(test du faisceau volé à la réponse de Dominic ).

R a une fonction intégrée dans le cadre de ses fonctions spéciales de mathématiques, y compris diverses formes de la fonction gamma.

Answer 4

Mathematica, 32 octets (30 caractères)

Sans aucune fonction Gammaintégrée associée, utilise la formule de Bubbler

Sum[#!/(-#2-x)^(#+1),{x,0,∞}]&

Essayez-le en ligne!

Mathematica, 29 octets

Sans PolyGamma[z]ouPolyGamma[n, z]

Log@Gamma@x~D~{x,#+1}/.x->#2&

Essayez-le en ligne!

Mathematica, 27 octets

Avec PolyGamma[z](c'est l'équivalent de la fonction digamma, ou \$\large\psi^0(z)\$)

PolyGamma@x~D~{x,#}/.x->#2&

Essayez-le en ligne!

Quelques programmes Mathematica qui n'utilisent pas le fichier intégré PolyGamma[n, z].

Answer 5

05AB1E , 16 15 octets

₄nÝ+I±mOI!IÉ·<P

-1 octet grâce à @ovs .

La première entrée est \$z\$, la deuxième entrée est \$m\$.

Essayez-le en ligne ou vérifiez tous les cas de test .

Explication:

Utilise le même algorithme que dans la réponse JavaScript de @Arnauld , alors assurez-vous de le voter.
Ou pour être plus précis, il utilise l'algorithme:

$$\psi^{(m)}(z)=(m\text{%}2\times2-1)\times m!\times\sum_{k=0}^{1000^2}{(z+k)^{\sim m}}$$

₄         # Push 1000
 n        # Square it to 1000000
  Ý       # Pop and push a list in the range [0,1000000]
   +      # Add the first (implicit) input-integer `z` to each value
    I     # Push the second input `m`
     ±    # Take it's bitwise-NOT: -m-1
      m   # Take each value to the power this `-m-1`
       O  # Sum all values in the list together
I!        # Push the second input `m` again, and take its factorial
IÉ        # Push the second input `m` again, and check if it's odd
          # (1 if truthy; 0 if falsey)
  ·       # Double that
   <      # And decrease it by 1
P         # And finally take the product of all three values on the stack
          # (after which it is output implicitly as result)

REMARQUE: s'il y a des inexactitudes très mineures dans les décimales, le ₄n( \$1\text{,}000\text{,}000\$) pourrait être remplacé par žm( \$9\text{,}876\text{,}543\text{,}210\$), même s'il serait trop lent de fonctionner sur TIO dans ce cas.

Answer 6

Octave / MATLAB, 4 octets

@psi

Essayez-le en ligne!

Answer 7

APL (Dyalog Unicode) 18,0, 20 octets

+/!⍤⊣÷(-(⍳!9)+⊢)*1+⊣

Essayez-le en ligne!

-2 octets grâce à Adám et ngn.

---

APL (Dyalog Unicode) , 22 octets

{+/(!⍺)÷(-⍵+⍳1e6)*1+⍺}

Essayez-le en ligne!

L'argument de gauche est \$m\$, l'argument droit est \$z\$.

Utilise une légère modification de la formule utilisée par d'autres réponses:

$$ \begin{aligned} \psi^{(m)}(z)&=(-1)^{m+1}m!\sum_{k=0}^{\infty}\dfrac{1}{(z+k)^{m+1}} \\ &\approx\sum_{k=0}^{10^6-1}\dfrac{m!}{(-z-k)^{m+1}} \end{aligned} $$

Comment ça fonctionne

{+/(!⍺)÷(-⍵+⍳1e6)*1+⍺}  ⍝ ⍺←m, ⍵←z
         -⍵+⍳1e6        ⍝ vector of -(z+0..999999)
        (       )*1+⍺   ⍝ raise each to the power of 1+m
 +/(!⍺)÷                ⍝ divide m! by each of above and sum them

Answer 8

R , 52 51 45 44 octets

Edit: -1 + 1 octets grâce à Giuseppe, qui a également souligné qu'il existe déjà une fonction R intégrée psigamma, qui résout la tâche pour seulement 8 octets

Edit2: ... et -6 octets de plus grâce à Robin Ryder

function(m,z)gamma(M<-m+1)*sum((-z:-1e4)^-M)

Essayez-le en ligne!

Utilise la même formule que la réponse d'Arnauld.
Les représentations de séries comme celle-ci sont très bien adaptées à R en tant que langage vectorisé nativement.

Remplacez les 1e4valeurs par des valeurs plus élevées (jusqu'à 9e9sans augmenter le nombre d'octets) pour une précision progressivement plus élevée et un temps d'exécution plus lent.

Answer 9

Python 3 , 38 octets

from scipy.special import*
f=polygamma

Essayez-le en ligne!

Answer 10

Java, 168 148 102 octets

(m,z)->{double p=1-m%2*2,f=0;long i=m;for(;i>0;)p*=i--;for(;i<1e7;)f-=p*Math.pow(z+i++,~m);return f;};

Explication

J'ai utilisé le même algorithme que dans la réponse JavaScript de @ Arnauld . Veuillez voter sa réponse.

Pour plus de commodité, voici la version formatée automatiquement:

(m, z) -> {
      double p = 1 - m % 2 * 2, f = 0;
      long i = m;
      for (; i > 0; ) p *= i--;
      for (; i < 1e7; ) f -= p * Math.pow(z + i++, ~m);
      return f;
    };

Code Java si typique: assez verbeux. Au moins ma version.

Edit: pourrait économiser 20 octets grâce à @user

Edit: économisé encore plus d'octets grâce à @ceilingcat

Essayez-le en ligne

Answer 11

Pyth , 24 octets

**^_1JhhQ*FhQsm^+deQ_JCG

Essayez-le en ligne! (le lien pointe vers un code légèrement différent qui additionne les 1e5termes au lieu de 1.56e62termes pour rendre le code exécutable et éviter les erreurs de débordement)

Explication

Utilise le même algorithme que dans la réponse JavaScript de @ Arnauld , alors assurez-vous de le voter.

**^_1JhhQ*FhQsm^+deQ_JCG
     JhhQ                  : Set J to first input + 1
  ^_1J                     : -1 ^ J
 *       *FhQ              : times factorial of first input
*            s             : times sum of
              m            : mapping
               ^+deQ_J     :   F(d): (d + (second input)) ^ -J
                      CG   : on range(1.56e62)

Answer 12

Scala , 68 octets

J'ai sauvé 3 octets et corrigé ma réponse grâce à Arnauld

Utilise l'algorithme de la réponse d'Arnauld

m=>z=>(0 to 1<<20 map(z+_ pow ~m)sum)*(m%2*2-1)*(1.0/:(1 to m))(_*_)

Essayez-le en ligne!