Renyi problem HW01 ICFP 2018.py

From Werner KRAUTH

Revision as of 00:01, 9 September 2018 Werner (Talk | contribs) (Renyi problem HW01 ICFP 2016.py moved to Renyi problem HW01 ICFP 2018.py) ← Previous diff		Current revision Werner (Talk | contribs)
Line 1:		Line 1:
-	This is the program Renyi.py that is useful for Homework 1 of the ICFP 2016 course in Statistical physics. Note that this program is written in Python 2.X.	+	This is the program Renyi.py that is useful for Homework 1 of the ICFP 2018 course in Statistical physics. Note that this program is written in Python 2.X. I hope to soon make a Python 3 program available.

---

Current revision

This is the program Renyi.py that is useful for Homework 1 of the ICFP 2018 course in Statistical physics. Note that this program is written in Python 2.X. I hope to soon make a Python 3 program available.

import random, math, scipy.misc, pylab

n = 5
x_val = []
renyi = []
int_renyi = []
for y in range(-200, 201):
    x = y / 10.0
    if abs(x) < n:
        kmax = int ((x + n)/ 2.0)
        renyi_sum = 0.0
        for k in range(0, kmax + 1):
            renyi_sum += (-1) ** k * scipy.misc.comb(n, k) * (n + x - 2 * k) ** (n - 1) /  2 **n / math.factorial(n - 1)
    else:
        renyi_sum = 0.0
    x_val.append(x)
    renyi.append(renyi_sum)
delx = x_val[1] - x_val[0]
for k in range(len(renyi)):
    int_renyi.append(sum(renyi[k:]) * delx)
pylab.ylim(0.0, 1.0)
pylab.plot(x_val, renyi, linewidth = 4, label='Renyi')
pylab.plot(x_val, int_renyi, linewidth = 4, label='Int Renyi')
pylab.legend(loc='upper right')
pylab.title("Renyi Formula for the sum of $n$ = " + str(n) + " indep.random ran(-1,1)")
pylab.xlabel('x')
pylab.ylabel('p(x)')
pylab.savefig('renyi.png')
pylab.show()