Code Monkey home page Code Monkey logo

orthopoly's Introduction

OrthoPoly

OrthoPoly is a Python class for generating orthogonal polynomials with respect to arbitrary probability density functions. The primary application of this script is in performing Arbitrary Polynomial Chaos Expansion (PCE) in uncertainty quantification studies.

Installation

Orthopoly is written in Python 3 and requires the following packages to be present in your system:

Since OrthoPoly is a simple class, no special installation steps are required. You can simply copy orthopoly.py to your working directory, and start using it.

Usage

OrthoPoly can be imported into your script as follows.

from orthopoly import OrthoPoly

To create an instance of OrthoPoly, you need to supply a probability density function (pdf). The code below defines the pdf of a random variable created by adding a uniform random variable and a normal random variable. This pdf is used to create an instance of OrthoPoly, and orthogonal polynomials are generated with respect to the supplied pdf.

def pdf(z, coeffs):
    mu, sigma, a, b = coeffs
    return 0.5/(b-a) * ( erf((z-a-mu)/sigma/sqrt(2)) - erf((z-b-mu)/sigma/sqrt(2)) )

pp = OrthoPoly(pdf, margs=[0, 1, -1, 1])
pp.gen_poly(5)

The function gen_poly() takes as an argument the largest order of the polynomial to be generated, and populates the variable pp.poly with the appropriate polynomials (which are numpy.polynomial.polynomial).

For numerical integration with respect to these polynomials, the quadrature points and weights can be obtained as follows.

points, weights = pp.get_quad_rule()

Numerical integration can also be performed directly using the quadrature() function. For instance, to compute the mean of the supplied pdf, one can do the following.

mean = pp.quadrature(lambda x: x)

For more examples, please look at the __main__ section of orthopoly.py, as well as the main.py script. main.py is a script which generates orthogonal polynomials for the pdf of a random variable which is a sum of a uniform and a normal random variable. It generates the polynomials, calculates the quadrature points and weights, and writes them to an output file.

References

To understand the mathematics behind OrthoPoly, please take a look at notes.pdf. For a more in-depth reading, consult [gautschi1982] and [golub1969].

[gautschi1982]Gautschi W. On Generating Orthogonal Polynomials. SIAM J Sci and Stat Comput 1982;3:289โ€“317. doi:10.1137/0903018.
[golub1969]Golub GH, Welsch JH. Calculation of Gauss Quadrature Rules. Mathematics of Computation 1969;23:221. doi:10.2307/2004418.

orthopoly's People

Contributors

j-jith avatar

Stargazers

 avatar  avatar  avatar  avatar  avatar

Watchers

 avatar  avatar

orthopoly's Issues

Severe instability - it is not possible to calculate more than ~10 polynomials

Example:

def pdf(x):
    return math.exp(-x)

pp = OrthoPoly(pdf, intlims=(0,1))
pp.gen_poly(n)

# a new class field 'betas' has been added to store the beta recurrence coefficients for demonstration
print(pp.betas)

Output:

[1.0,
0.07932640579220766,
0.06711669572587618,
0.06433876970545459,
0.06350596398737984,
0.06313659912953457,
0.06293951687236084,
0.06282180581953442,
0.0627458612867561,
0.06269478868364445,
0.06283407650853606,
0.06089373421606997,
0.21595648578197688,
0.669882586586399,
-1.3444323914957315,
4.962392290085074,
0.4906144952293512,
0.20512208494961423,
-686.035257250946,
15925.79931665891,
56662.12826844703]

The last elements demonstrate instability.

Recommend Projects

  • React photo React

    A declarative, efficient, and flexible JavaScript library for building user interfaces.

  • Vue.js photo Vue.js

    ๐Ÿ–– Vue.js is a progressive, incrementally-adoptable JavaScript framework for building UI on the web.

  • Typescript photo Typescript

    TypeScript is a superset of JavaScript that compiles to clean JavaScript output.

  • TensorFlow photo TensorFlow

    An Open Source Machine Learning Framework for Everyone

  • Django photo Django

    The Web framework for perfectionists with deadlines.

  • D3 photo D3

    Bring data to life with SVG, Canvas and HTML. ๐Ÿ“Š๐Ÿ“ˆ๐ŸŽ‰

Recommend Topics

  • javascript

    JavaScript (JS) is a lightweight interpreted programming language with first-class functions.

  • web

    Some thing interesting about web. New door for the world.

  • server

    A server is a program made to process requests and deliver data to clients.

  • Machine learning

    Machine learning is a way of modeling and interpreting data that allows a piece of software to respond intelligently.

  • Game

    Some thing interesting about game, make everyone happy.

Recommend Org

  • Facebook photo Facebook

    We are working to build community through open source technology. NB: members must have two-factor auth.

  • Microsoft photo Microsoft

    Open source projects and samples from Microsoft.

  • Google photo Google

    Google โค๏ธ Open Source for everyone.

  • D3 photo D3

    Data-Driven Documents codes.