StackOverflow2013

Note that there are some explanatory texts on larger screens.

plurals

PO
text
Body
copied!<p>If I understand your question correctly, you want to calculate the integral of f(x) * Ln(x) where f(x) is a piecewise function you're defining with a python function. I'm assuming you're not specifically interested in this particular step function.</p> <p>You can get the values of the Legendre polynomials using <a href="http://docs.scipy.org/doc/numpy/reference/generated/numpy.polynomial.legendre.legval.html#numpy.polynomial.legendre.legval" rel="nofollow noreferrer">legval</a> with the identity matrix for the coefficient argument.</p> <pre><code>import numpy as np import matplotlib x = np.linspace(-1, 1, 201) L = np.polynomial.legendre.legval(x, np.identity(50)) plt.plot(x, L.T) </code></pre> <p><img src="https://i.stack.imgur.com/H8KpJ.png" alt="enter image description here"></p> <p>You can then perform the integral with quadrature. Using gauss-legendre quadrature might be more efficient since the integral of a legendre polynomial will be exact for Ln(x) where n is less than the quadrature size.</p> <pre><code>import numpy as np from numpy.polynomial.legendre import leggauss, legval def f(x): if 0 <= x <= 1: return 1 if -1 <= x <= 0: return -1 # of course you could write a vectorized version of # this particular f(x), but I assume you have a more # general piecewise function f = np.vectorize(f) deg = 100 x, w = leggauss(deg) # len(x) == 100 L = np.polynomial.legendre.legval(x, np.identity(deg)) # Sum L(xi)*f(xi)*wi integral = (L*(f(x)*w)[None,:]).sum(axis=1) c = (np.arange(1,51) + 0.5) * integral[1:51] x_fine = np.linspace(-1, 1, 2001) # 2001 points Lfine = np.polynomial.legendre.legval(x_fine, np.identity(51)) # sum_1_50 of c(n) * Ln(x_fine) cLn_sum = (c[:,None] * Lfine[1:51,:]).sum(axis=0) </code></pre> <p>c = 1.5, 0, -8.75e-1, 0, ... which I think is the result you're looking for.</p>

Querying!

Guidance

An individual column

Larger individual text columns get their own page to allow for proper reading.

SQuiL has stopped working due to an internal error.

If you are curious you may find further information in the browser console, which is accessible through the devtools (F12).

Reload