<SYSTEM_TASK:>
Solve the following problem using Python, implementing the functions described below, one line at a time
<END_TASK>
<USER_TASK:>
Description:
 def qnwequi(n, a, b, kind="N", equidist_pp=None, random_state=None):
 """ Generates equidistributed sequences with property that averages value of integrable function evaluated over the sequence converges to the integral as n goes to infinity. Parameters n : int Number of sequence points a : scalar or array_like(float) A length-d iterable of lower endpoints. If a scalar is given, that constant is repeated d times, where d is the number of dimensions b : scalar or array_like(float) A length-d iterable of upper endpoints. If a scalar is given, that constant is repeated d times, where d is the number of dimensions kind : string, optional(default="N") One of the following: - N - Neiderreiter (default) - W - Weyl - H - Haber - R - pseudo Random equidist_pp : array_like, optional(default=None) TODO: I don't know what this does random_state : int or np.random.RandomState, optional Random seed (integer) or np.random.RandomState instance to set the initial state of the random number generator for reproducibility. If None, a randomly initialized RandomState is used. Returns ------- nodes : np.ndarray(dtype=float) Quadrature nodes weights : np.ndarray(dtype=float) Weights for quadrature nodes Notes ----- Based of original function ``qnwequi`` in CompEcon toolbox by Miranda and Fackler References Miranda, Mario J, and Paul L Fackler. Applied Computational Economics and Finance, MIT Press, 2002. """

    random_state = check_random_state(random_state)

    if equidist_pp is None:
        import sympy as sym
        equidist_pp = np.sqrt(np.array(list(sym.primerange(0, 7920))))

    n, a, b = list(map(np.atleast_1d, list(map(np.asarray, [n, a, b]))))

    d = max(list(map(len, [n, a, b])))
    n = np.prod(n)

    if a.size == 1:
        a = np.repeat(a, d)

    if b.size == 1:
        b = np.repeat(b, d)

    i = np.arange(1, n + 1)

    if kind.upper() == "N":  # Neiderreiter
        j = 2.0 ** (np.arange(1, d+1) / (d+1))
        nodes = np.outer(i, j)
        nodes = (nodes - fix(nodes)).squeeze()
    elif kind.upper() == "W":  # Weyl
        j = equidist_pp[:d]
        nodes = np.outer(i, j)
        nodes = (nodes - fix(nodes)).squeeze()
    elif kind.upper() == "H":  # Haber
        j = equidist_pp[:d]
        nodes = np.outer(i * (i+1) / 2, j)
        nodes = (nodes - fix(nodes)).squeeze()
    elif kind.upper() == "R":  # pseudo-random
        nodes = random_state.rand(n, d).squeeze()
    else:
        raise ValueError("Unknown sequence requested")

    # compute nodes and weights
    r = b - a
    nodes = a + nodes * r
    weights = (np.prod(r) / n) * np.ones(n)

    return nodes, weights