Generating numbers in log-normal distribution not giving right parameters in return?

So I am trying to implement a function to generate random number following log-normal distribution myself. (Despite I am testing it in numpy, I need it for somewhere else where out of the box function is unavailable).

Using Box-Muller transform to generate Normally distributed random number and taking exponential of it (suggested in this question), and referring to Wikipedia page I come up with this code:

import matplotlib.pyplot as plt
import numpy as np

rng = np.random.default_rng()

def f():
    mean = 100
    sigma = 50

    u = np.log((mean * mean) / np.sqrt(mean * mean + sigma * sigma))
    d = np.log(1 + (sigma * sigma) / (mean * mean))

    u1 = 0
    while u1 <= 0.0000001:
        u1 = rng.random()
    u2 = rng.random()

    mag = d * np.sqrt(-2 * np.log(u1))
    z0 = mag * np.cos(2 * np.pi * u2) + u
    z1 = mag * np.sin(2 * np.pi * u2) + u

    return (np.exp(z0), np.exp(z1))

iteration = 100000

data = []
for i in range(0, iteration):
    x = f()
    data.append(x[0])
    data.append(x[1])

fig, ax = plt.subplots()
ax.hist(data, bins=100) 

print('mean: ' + str(np.mean(data)))
print('median: ' + str(np.median(data)))
print('stdev: ' + str(np.std(data)))

However, the generated values does not have the right mean and stdev, even though the shape seems to be vaguely correct. Does I miss anything? Or is it just some floating point shenanigans?