Simulace Cauchyho rozdělení metodou inverzní funkce

$r\in U(0,1)$
transformace: $x = \text{tg}\left (\pi \left (r-\frac{1}{2}\right )\right )$

import numpy as np
import matplotlib.pyplot as plt

def Cauchy(x):
    return 1/np.pi*1/(1+x**2)

rng=np.random.Generator(np.random.PCG64())

N=1000000 #number of simulations
r=rng.random(N)
x=np.tan(np.pi*(r-0.5))

xmin=-100 #range for histrogram - min
xmax=100  #range for histrogram - max
nbins=5000 #number of bins
plt.hist(x,range=[xmin,xmax],bins=nbins,density=True)
x_plot=np.linspace(xmin,xmax,nbins)
plt.plot(x_plot,Cauchy(x_plot),c='red')
plt.xlabel('x')
plt.ylabel('f(x)')
plt.xlim(-5,5)
plt.show()