def legendre( order, lower=-1, upper=1, physicist=False, normed=False, retall=False, ): """ Examples: >>> polynomials, norms = chaospy.expansion.legendre(3, retall=True) >>> polynomials.round(5) polynomial([1.0, q0, q0**2-0.33333, q0**3-0.6*q0]) >>> norms array([1. , 0.33333333, 0.08888889, 0.02285714]) >>> chaospy.expansion.legendre(3, physicist=True).round(3) polynomial([1.0, 1.5*q0, 2.5*q0**2-0.556, 4.375*q0**3-1.672*q0]) >>> chaospy.expansion.legendre(3, lower=0, upper=1, normed=True).round(3) polynomial([1.0, 3.464*q0-1.732, 13.416*q0**2-13.416*q0+2.236, 52.915*q0**3-79.373*q0**2+31.749*q0-2.646]) """ multiplier = 1.0 if physicist: multiplier = numpy.arange(1, order + 1) multiplier = (2 * multiplier + 1) / (multiplier + 1) _, [polynomials], [norms] = chaospy.recurrence.analytical_stieltjes( order, chaospy.Uniform(lower, upper), multiplier=multiplier) if normed: polynomials = chaospy.true_divide(polynomials, numpy.sqrt(norms)) norms[:] = 1.0 return (polynomials, norms) if retall else polynomials
def laguerre( order, alpha=0.0, physicist=False, normed=False, retall=False, ): """ Examples: >>> polynomials, norms = chaospy.expansion.laguerre(3, retall=True) >>> polynomials polynomial([1.0, q0-1.0, q0**2-4.0*q0+2.0, q0**3-9.0*q0**2+18.0*q0-6.0]) >>> norms array([ 1., 1., 4., 36.]) >>> chaospy.expansion.laguerre(3, physicist=True).round(5) polynomial([1.0, -q0+1.0, 0.5*q0**2-2.0*q0+2.0, -0.16667*q0**3+1.5*q0**2-5.33333*q0+4.66667]) >>> chaospy.expansion.laguerre(3, alpha=2, normed=True).round(3) polynomial([1.0, 0.577*q0-1.732, 0.204*q0**2-1.633*q0+2.449, 0.053*q0**3-0.791*q0**2+3.162*q0-3.162]) """ multiplier = -1.0 / numpy.arange(1, order + 1) if physicist else 1.0 _, [polynomials], [norms] = chaospy.recurrence.analytical_stieltjes( order, chaospy.Gamma(alpha + 1), multiplier=multiplier) if normed: polynomials = chaospy.true_divide(polynomials, numpy.sqrt(norms)) norms[:] = 1.0 return (polynomials, norms) if retall else polynomials
def gegenbauer( order, alpha, lower=-1, upper=1, physicist=False, normed=False, retall=False, ): """ Gegenbauer polynomials. Args: order (int): The polynomial order. alpha (float): Gegenbauer shape parameter. lower (float): Lower bound for the integration interval. upper (float): Upper bound for the integration interval. physicist (bool): Use physicist weights instead of probabilist. Examples: >>> polynomials, norms = chaospy.expansion.gegenbauer(4, 1, retall=True) >>> polynomials polynomial([1.0, q0, q0**2-0.25, q0**3-0.5*q0, q0**4-0.75*q0**2+0.0625]) >>> norms array([1. , 0.25 , 0.0625 , 0.015625 , 0.00390625]) >>> chaospy.expansion.gegenbauer(3, 1, physicist=True) polynomial([1.0, 2.0*q0, 4.0*q0**2-0.5, 8.0*q0**3-2.0*q0]) >>> chaospy.expansion.gegenbauer(3, 1, lower=0.5, upper=1.5, normed=True).round(3) polynomial([1.0, 4.0*q0-4.0, 16.0*q0**2-32.0*q0+15.0, 64.0*q0**3-192.0*q0**2+184.0*q0-56.0]) """ multiplier = 1 if physicist: multiplier = numpy.arange(1, order + 1) multiplier = 2 * (multiplier + alpha - 1) / multiplier _, [polynomials], [norms] = chaospy.recurrence.analytical_stieltjes( order, chaospy.Beta(alpha + 0.5, alpha + 0.5, lower, upper), multiplier=multiplier, ) if normed: polynomials = chaospy.true_divide(polynomials, numpy.sqrt(norms)) norms[:] = 1.0 return (polynomials, norms) if retall else polynomials
def hermite( order, mu=0.0, sigma=1.0, physicist=False, normed=False, retall=False, ): """ Hermite orthogonal polynomial expansion. Args: order (int): The quadrature order. mu (float): Non-centrality parameter. sigma (float): Scale parameter. physicist (bool): Use physicist weights instead of probabilist variant. normed (bool): If True orthonormal polynomials will be used. retall (bool): If true return numerical stabilized norms as well. Roughly the same as ``cp.E(orth**2, dist)``. Returns: (numpoly.ndpoly, numpy.ndarray): Hermite polynomial expansion. Norms of the orthogonal expansion on the form ``E(orth**2, dist)``. Examples: >>> polynomials, norms = chaospy.expansion.hermite(4, retall=True) >>> polynomials polynomial([1.0, q0, q0**2-1.0, q0**3-3.0*q0, q0**4-6.0*q0**2+3.0]) >>> norms array([ 1., 1., 2., 6., 24.]) >>> chaospy.expansion.hermite(3, physicist=True) polynomial([1.0, 2.0*q0, 4.0*q0**2-2.0, 8.0*q0**3-12.0*q0]) >>> chaospy.expansion.hermite(3, sigma=2.5, normed=True).round(3) polynomial([1.0, 0.4*q0, 0.113*q0**2-0.707, 0.026*q0**3-0.49*q0]) """ multiplier = 2 if physicist else 1 _, [polynomials], [norms] = chaospy.recurrence.analytical_stieltjes( order, chaospy.Normal(mu, sigma), multiplier=multiplier) if normed: polynomials = chaospy.true_divide(polynomials, numpy.sqrt(norms)) norms[:] = 1.0 return (polynomials, norms) if retall else polynomials
def chebyshev_2( order, lower=-1, upper=1, physicist=False, normed=False, retall=False, ): """ Chebyshev polynomials of the second kind. Args: order (int): The quadrature order. lower (float): Lower bound for the integration interval. upper (float): Upper bound for the integration interval. physicist (bool): Use physicist weights instead of probabilist. Returns: (numpoly.ndpoly, numpy.ndarray): Chebyshev polynomial expansion. Norms of the orthogonal expansion on the form ``E(orth**2, dist)``. Examples: >>> polynomials, norms = chaospy.expansion.chebyshev_2(4, retall=True) >>> polynomials polynomial([1.0, q0, q0**2-0.25, q0**3-0.5*q0, q0**4-0.75*q0**2+0.0625]) >>> norms array([1. , 0.25 , 0.0625 , 0.015625 , 0.00390625]) >>> chaospy.expansion.chebyshev_2(3, physicist=True) polynomial([1.0, 2.0*q0, 4.0*q0**2-0.5, 8.0*q0**3-2.0*q0]) >>> chaospy.expansion.chebyshev_2(3, lower=0.5, upper=1.5, normed=True).round(3) polynomial([1.0, 4.0*q0-4.0, 16.0*q0**2-32.0*q0+15.0, 64.0*q0**3-192.0*q0**2+184.0*q0-56.0]) """ multiplier = 2 if physicist else 1 _, [polynomials], [norms] = chaospy.recurrence.analytical_stieltjes( order, chaospy.Beta(1.5, 1.5, lower, upper), multiplier=multiplier ) if normed: polynomials = chaospy.true_divide(polynomials, numpy.sqrt(norms)) norms[:] = 1.0 return (polynomials, norms) if retall else polynomials
def chebyshev_1( order, lower=-1, upper=1, physicist=False, normed=False, retall=False, ): """ Chebyshev polynomials of the first kind. Args: order (int): The polynomial order. lower (float): Lower bound for the integration interval. upper (float): Upper bound for the integration interval. physicist (bool): Use physicist weights instead of probabilist. Returns: (numpoly.ndpoly, numpy.ndarray): Chebyshev polynomial expansion. Norms of the orthogonal expansion on the form ``E(orth**2, dist)``. Examples: >>> polynomials, norms = chaospy.expansion.chebyshev_1(4, retall=True) >>> polynomials polynomial([1.0, q0, q0**2-0.5, q0**3-0.75*q0, q0**4-q0**2+0.125]) >>> norms array([1. , 0.5 , 0.125 , 0.03125 , 0.0078125]) >>> chaospy.expansion.chebyshev_1(3, physicist=True) polynomial([1.0, q0, 2.0*q0**2-1.0, 4.0*q0**3-2.5*q0]) >>> chaospy.expansion.chebyshev_1(3, lower=0.5, upper=1.5, normed=True).round(3) polynomial([1.0, 2.828*q0-2.828, 11.314*q0**2-22.627*q0+9.899, 45.255*q0**3-135.765*q0**2+127.279*q0-36.77]) """ multiplier = 1 + numpy.arange(order).astype(bool) if physicist else 1 _, [polynomials], [norms] = chaospy.recurrence.analytical_stieltjes( order, chaospy.Beta(0.5, 0.5, lower, upper), multiplier=multiplier ) if normed: polynomials = chaospy.true_divide(polynomials, numpy.sqrt(norms)) norms[:] = 1.0 return (polynomials, norms) if retall else polynomials
def jacobi( order, alpha, beta, lower=-1, upper=1, physicist=False, normed=False, retall=False, ): """ Jacobi polynomial expansion. Examples: >>> polynomials, norms = chaospy.expansion.jacobi(4, 0.5, 0.5, retall=True) >>> polynomials polynomial([1.0, q0, q0**2-0.5, q0**3-0.75*q0, q0**4-q0**2+0.125]) >>> norms array([1. , 0.5 , 0.125 , 0.03125 , 0.0078125]) >>> chaospy.expansion.jacobi(3, 0.5, 0.5, physicist=True).round(4) polynomial([1.0, 1.5*q0, 2.5*q0**2-0.8333, 4.375*q0**3-2.1146*q0]) >>> chaospy.expansion.jacobi(3, 1.5, 0.5, normed=True) polynomial([1.0, 2.0*q0, 4.0*q0**2-1.0, 8.0*q0**3-4.0*q0]) """ multiplier = 1 if physicist: multiplier = numpy.arange(1, order + 1) multiplier = ((2 * multiplier + alpha + beta - 1) * (2 * multiplier + alpha + beta) / (2 * multiplier * (multiplier + alpha + beta))) _, [polynomials], [norms] = chaospy.recurrence.analytical_stieltjes( order, chaospy.Beta(alpha, beta, lower=lower, upper=upper), multiplier=multiplier, ) if normed: polynomials = chaospy.true_divide(polynomials, numpy.sqrt(norms)) norms[:] = 1.0 return (polynomials, norms) if retall else polynomials