def compute_b_mp(sigma, q, lmbd, verbose=False):
    lmbd_int = int(math.ceil(lmbd))
    if (lmbd_int == 0):
        return 1.0
    (mu0, _, mu) = distributions_mp(sigma, q)
    b_lambda_fn = (lambda z: (mu0(z) * ((mu0(z) / mu(z)) ** lmbd_int)))
    b_lambda = integral_inf_mp(b_lambda_fn)
    m = ((sigma ** 2) * (mp.log(((2 - q) / (1 - q))) + (1 / (2 * (sigma ** 2)))))
    b_fn = (lambda z: (((mu0(z) / mu(z)) ** lmbd_int) - ((mu((- z)) / mu0(z)) ** lmbd_int)))
    if verbose:
        print('M =', m)
        print('f(-M) = {} f(M) = {}'.format(b_fn((- m)), b_fn(m)))
        assert ((b_fn((- m)) < 0) and (b_fn(m) < 0))
    b_lambda_int1_fn = (lambda z: (mu0(z) * ((mu0(z) / mu(z)) ** lmbd_int)))
    b_lambda_int2_fn = (lambda z: (mu0(z) * ((mu(z) / mu0(z)) ** lmbd_int)))
    b_int1 = integral_bounded_mp(b_lambda_int1_fn, (- m), m)
    b_int2 = integral_bounded_mp(b_lambda_int2_fn, (- m), m)
    a_lambda_m1 = compute_a_mp(sigma, q, (lmbd - 1))
    b_bound = ((a_lambda_m1 + b_int1) - b_int2)
    if verbose:
        print('B by numerical integration', b_lambda)
        print('B must be no more than    ', b_bound)
    assert (b_lambda < (b_bound + 1e-05))
    return _to_np_float64(b_lambda)
def compute_b_mp(sigma, q, lmbd, verbose=False):
    lmbd_int = int(math.ceil(lmbd))
    if (lmbd_int == 0):
        return 1.0
    (mu0, _, mu) = distributions_mp(sigma, q)
    b_lambda_fn = (lambda z: (mu0(z) * ((mu0(z) / mu(z)) ** lmbd_int)))
    b_lambda = integral_inf_mp(b_lambda_fn)
    m = ((sigma ** 2) * (mp.log(((2 - q) / (1 - q))) + (1 / (2 * (sigma ** 2)))))
    b_fn = (lambda z: (((mu0(z) / mu(z)) ** lmbd_int) - ((mu((- z)) / mu0(z)) ** lmbd_int)))
    if verbose:
        print('M =', m)
        print('f(-M) = {} f(M) = {}'.format(b_fn((- m)), b_fn(m)))
        assert ((b_fn((- m)) < 0) and (b_fn(m) < 0))
    b_lambda_int1_fn = (lambda z: (mu0(z) * ((mu0(z) / mu(z)) ** lmbd_int)))
    b_lambda_int2_fn = (lambda z: (mu0(z) * ((mu(z) / mu0(z)) ** lmbd_int)))
    b_int1 = integral_bounded_mp(b_lambda_int1_fn, (- m), m)
    b_int2 = integral_bounded_mp(b_lambda_int2_fn, (- m), m)
    a_lambda_m1 = compute_a_mp(sigma, q, (lmbd - 1))
    b_bound = ((a_lambda_m1 + b_int1) - b_int2)
    if verbose:
        print('B by numerical integration', b_lambda)
        print('B must be no more than    ', b_bound)
    assert (b_lambda < (b_bound + 1e-05))
    return _to_np_float64(b_lambda)