def qbfs_recurrence_P(n, x, Pnm1=None, Pnm2=None, recursion_coef=None):
    """P(m+1) from oe-18-19-19700 eq. (2.6).

    Parameters
    ----------
    n : `int`
        polynomial order
    x : `numpy.ndarray`
        x values, notionally in / orthogonal over [0, 1], to evaluate at
    Pnm1 : `numpy.ndarray`, optional
        the value of this function for argument n - 1
    Pnm2 : `numpy.ndarray`, optional
        the value of this function for argument n - 2
    recursion_coef : `numpy.ndarray`, optional
        the coefficient to apply, if recursion_coef = C: evaluates C * Pnm1 - Pnm2

    Returns
    -------
    `numpy.ndarray`
        the value of the auxiliary P polynomial for given order n and point(s) x

    """
    if n == 0:
        return 2
    elif n == 1:
        return 6 - 8 * x
    else:
        if Pnm2 is None:
            Pnm2 = qbfs_recurrence_P(n - 2, x)
        if Pnm1 is None:
            Pnm1 = qbfs_recurrence_P(n - 1, x, Pnm1=Pnm2)

        if recursion_coef is None:
            recursion_coef = 2 - 4 * x

        return recursion_coef * Pnm1 - Pnm2