Implementa la funzione Polygamma

Sep 01 2020

La funzione Polygamma di ordine \$m\$, \$\psi^{(m)}(z)\$, è il \$(m + 1)\$esima derivata del logaritmo della funzione gamma , che è anche la \$m\$esima derivata della funzione digamma . Il tuo compito è prendere un numero intero \$m\$e un numero reale positivo \$z\$e output \$\psi^{(m)}(z)\$

Definizioni

Per chi non ha familiarità con le funzioni di cui sopra (Gamma, Digamma e Polygamma), ecco alcune definizioni diverse per ciascuna:

\$\Gamma(z)\$

  • La funzione gamma è un'estensione del fattoriale ( \$x! = 1\cdot2\cdot3\cdots(x-1)\cdot(x)\$) a numeri reali
  • \$\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 funzione digamma è la derivata logaritmica della funzione 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 funzione polygamma di ordine \$m\$è il \$m\$esima derivata della funzione 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}}\$

Compito

Devi prendere due input, un numero naturale \$m\$e un numero reale positivo \$z\$e output \$\psi^{(m)}(z)\$. Gli input e gli output rientreranno sempre nei limiti numerici della tua lingua, ma il tuo algoritmo deve funzionare teoricamente per tutti gli input .

Poiché l'output di solito sarà un numero reale, anziché un intero, l'output dovrebbe essere corretto con almeno 10 cifre significative. Gli zeri finali possono essere omessi per i valori esatti. Ad esempio, se l'output è un numero intero, i decimali finali 0non sono obbligatori, ma sono consentiti se lo desideri.

Questo è il codice del golf, quindi vince il codice più breve in byte.

Casi test

I risultati possono differire a causa di imprecisioni in virgola mobile, la libreria scipy di Python è stata utilizzata per generare i valori. I valori sono arrotondati a 15d.p., salvo diversa indicazione.

 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

Questa tabella ha i valori di \$\psi^{(m)}(z)\$per \$0 \le m \le 9\$e \$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  |
+---+------------------------+---------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+--------------------+

Risposte

16 Arnauld Sep 01 2020 at 14:17

JavaScript (ES7),  68 66 61  59 byte

Aspetta (m)(z).

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

Provalo online!

Questo si basa sulla seguente rappresentazione in serie ( da Wikipedia ):

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

Commentato

(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
12 user Sep 01 2020 at 16:21

Mathematica, 9 byte

PolyGamma

Ovviamente Mathematica ha un builtin per questo

8 Giuseppe Sep 01 2020 at 17:01

R , 8 byte

psigamma

Provalo online!

Prende input z,m(test harness rubato dalla risposta di Dominic ).

R ha un builtin come parte delle sue funzioni speciali di matematica incluse varie forme della funzione gamma.

8 ZippyMagician Sep 01 2020 at 19:56

Mathematica, 32 byte (30 caratteri)

Senza alcun Gammabuiltin correlato, utilizza la formula di Bubbler

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

Provalo online!

Mathematica, 29 byte

Senza PolyGamma[z]oPolyGamma[n, z]

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

Provalo online!

Mathematica, 27 byte

Con PolyGamma[z](questo è l'equivalente della funzione digamma, o \$\large\psi^0(z)\$)

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

Provalo online!

Alcuni programmi Mathematica che non utilizzano il builtin PolyGamma[n, z].

7 KevinCruijssen Sep 01 2020 at 15:04

05AB1E , 16 15 byte

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

-1 byte grazie a @ovs .

Il primo input è \$z\$, il secondo input è \$m\$.

Provalo online o verifica tutti i casi di test .

Spiegazione:

Utilizza lo stesso algoritmo della risposta JavaScript di @Arnauld , quindi assicurati di votarlo.
O per essere più precisi, utilizza l'algoritmo:

$$\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)

NOTA: se sono presenti inesattezze molto minori nei decimali, il ₄n( \$1\text{,}000\text{,}000\$) potrebbe essere sostituito con žm( \$9\text{,}876\text{,}543\text{,}210\$), anche se in tal caso sarebbe troppo lento eseguire su TIO.

5 LuisMendo Sep 01 2020 at 17:25

Ottava / MATLAB, 4 byte

@psi

Provalo online!

5 Bubbler Sep 02 2020 at 01:13

APL (Dyalog Unicode) 18.0, 20 byte

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

Provalo online!

-2 byte grazie ad Adám e ngn.


APL (Dyalog Unicode) , 22 byte

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

Provalo online!

L'argomento a sinistra è \$m\$, l'argomento giusto è \$z\$.

Utilizza una leggera modifica della formula utilizzata da altre risposte:

$$ \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} $$

Come funziona

{+/(!⍺)÷(-⍵+⍳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
4 DominicvanEssen Sep 01 2020 at 16:03

R , 52 51 45 44 byte

Modifica: -1 + 1 byte grazie a Giuseppe, che ha anche sottolineato che esiste già una funzione R incorporata psigamma, che risolve il compito per soli 8 byte

Edit2: ... e altri -6 byte grazie a Robin Ryder

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

Provalo online!

Utilizza la stessa formula della risposta di Arnauld.
Le rappresentazioni in serie come questa sono molto adatte a R come linguaggio a vettorializzazione nativa.

Modificare il 1e4valore in valori più alti (fino a 9e9senza aumentare il numero di byte) per una precisione progressivamente maggiore e un tempo di esecuzione più lento.

2 Noodle9 Sep 01 2020 at 22:59

Python 3 , 38 byte

from scipy.special import*
f=polygamma

Provalo online!

2 jmizv Sep 03 2020 at 10:27

Java, 168 148102 byte

(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;};

Spiegazione

Ho usato lo stesso algoritmo della risposta JavaScript di @ Arnauld . Si prega di votare la sua risposta.

Per comodità, ecco la versione formattata automaticamente:

(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;
    };

Così tipico codice Java: piuttosto prolisso. Almeno la mia versione.

Modifica: potrebbe salvare 20 byte grazie a @user

Modifica: salvato ancora più byte grazie a @ceilingcat

Provalo online

1 Mukundan314 Sep 01 2020 at 15:56

Pyth , 24 byte

**^_1JhhQ*FhQsm^+deQ_JCG

Provalo online! (il collegamento punta a un codice leggermente diverso che somma i 1e5termini anziché i 1.56e62termini per rendere il codice eseguibile ed evitare errori di overflow)

Spiegazione

Utilizza lo stesso algoritmo della risposta JavaScript di @ Arnauld , quindi assicurati di dargli un voto positivo.

**^_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)
1 user Sep 01 2020 at 16:03

Scala , 68 byte

Ho salvato 3 byte e ho corretto la mia risposta grazie ad Arnauld

Utilizza l'algoritmo dalla risposta di Arnauld

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

Provalo online!