2
0
mirror of https://github.com/boostorg/math.git synced 2026-02-21 15:12:28 +00:00

[CI SKIP] Use JM version of hypergeometric.gbk

This commit is contained in:
pabristow
2019-08-09 16:13:10 +01:00
5 changed files with 149 additions and 2 deletions

View File

@@ -30,6 +30,8 @@ LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
Real operator()(Real x) const;
Real prime(Real x) const;
Real double_prime(Real x) const;
};
}} // namespaces
@@ -37,7 +39,7 @@ LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
[heading Cubic B-Spline Interpolation]
The cubic B-spline class provided by boost allows fast and accurate interpolation of a function which is known at equally spaced points.
The cubic B-spline class provided by Boost allows fast and accurate interpolation of a function which is known at equally spaced points.
The cubic B-spline interpolation is numerically stable as it uses compactly supported basis functions constructed via iterative convolution.
This is to be contrasted to traditional cubic spline interpolation which is ill-conditioned as the global support of cubic polynomials causes small changes far from the evaluation point to exert a large influence on the calculated value.
@@ -81,6 +83,15 @@ and to evaluate the derivative of the interpolant we use
Be aware that the accuracy guarantees on the derivative of the spline are an order lower than the guarantees on the original function,
see [@http://www.springer.com/us/book/9780387984087 Numerical Analysis, Graduate Texts in Mathematics, 181, Rainer Kress] for details.
The last interesting member is the second derivative, evaluated via
double ypp = spline.double_prime(x);
The basis functions of the spline are cubic polynomials, so the second derivative is simply linear interpolation.
But the interpolation is not constrained by data (you don't pass in the second derivatives), and hence is less accurate than would be expected from assuming it is a linear interpolation.
The problem is especially pronounced at the boundaries, where the second derivative is totally unconstrained.
Use the second derivative of the cubic B-spline interpolator only in desperation.
Finally, note that this is an interpolator, not an extrapolator.
Therefore, you should strenuously avoid evaluating the spline outside the endpoints.
However, it is not an error if you do, as often you cannot control where (say) an ODE stepper will evaluate your function.

View File

@@ -168,6 +168,7 @@ Otherwise we use the defining series.
The function `hypergeometric_1F1(a, b, z)` returns the non-singular solution to
[@https://en.wikipedia.org/wiki/Confluent_hypergeometric_function Kummer's equation]
<<<<<<< HEAD
[/\large $$z \frac{d^2 w}{d z^2} + (b-z) \frac{dw}{dz} - aw = 0 $$]
[equation hypergeometric_1f1_2]
@@ -176,12 +177,21 @@ The function `hypergeometric_1F1(a, b, z)` returns the non-singular solution to
which for |/z/| < 1 has the hypergeometric series expansion
[equation hypergeometric_1f1_1]
=======
[:[/\large $$z \frac{d^2 w}{d z^2} + (b-z) \frac{dw}{dz} - aw = 0 $$]
[$../equations/hypergeometric_1f1_2.svg]]
which for |/z/| < 1 has the hypergeometric series expansion
[:[$../equations/hypergeometric_1f1_1.svg]]
>>>>>>> 75789acd778461428be66b3bbc3ada25e33a042c
where (/a/)[sub /n/] denotes rising factorial.
This function has the same definition as Mathematica's `Hypergeometric1F1[a, b, z]` and Maple's `KummerM(a, b, z)`.
The "regularized" versions of the function return:
<<<<<<< HEAD
[/ \Large $$ \textbf{M}(a, b; z) = \frac{{_1F_1}(a, b; z)}{\Gamma(b)} = \sum_{n=0}^{\infty} \frac{(a)_n z^n}{\Gamma(b+n) n!} $$]
[equation hypergeometric_1f1_17]
@@ -189,6 +199,15 @@ The "log" versions of the function return:
[/ \Large $$ \ln(|_1F_1(a, b, z)|) $$ ]
[equation hypergeometric_1f1_18]
=======
[:[/ \Large $$ \textbf{M}(a, b; z) = \frac{{_1F_1}(a, b; z)}{\Gamma(b)} = \sum_{n=0}^{\infty} \frac{(a)_n z^n}{\Gamma(b+n) n!} $$]
[$../equations/hypergeometric_1f1_17.svg]]
The "log" versions of the function return:
[:[/ \Large $$ \ln(|_1F_1(a, b, z)|) $$ ]
[$../equations/hypergeometric_1f1_18.svg]]
>>>>>>> 75789acd778461428be66b3bbc3ada25e33a042c
When the optional `sign` argument is supplied it is set on exit to either +1 or -1 depending on the sign of [sub 1]/F/[sub 1](/a/, /b/, /z/).
@@ -277,12 +296,21 @@ readers will have to refer to the code for the transitions between methods, as w
In almost all cases where /z < 0/ we use [@https://en.wikipedia.org/wiki/Confluent_hypergeometric_function Kummer's relation]
to make /z/ positive:
<<<<<<< HEAD
[/\large $$_1F_1(a, b; -z) = e^{-z} {_1F_1}(a, b; z)$$]
[equation hypergeometric_1f1_12]
The main series representation
[equation hypergeometric_1f1_1]
=======
[:[/\large $$_1F_1(a, b; -z) = e^{-z} {_1F_1}(a, b; z)$$]
[$../equations/hypergeometric_1f1_12.svg]]
The main series representation
[:[$../equations/hypergeometric_1f1_1.svg]]
>>>>>>> 75789acd778461428be66b3bbc3ada25e33a042c
is used only when
@@ -294,13 +322,19 @@ to life again as /b/ crosses the origin.
A&S 13.3.6 gives
<<<<<<< HEAD
[/\large $$_1F_1(a, b; z) = e^{ \frac{z}{2} } \Gamma(b-a- \frac{1}{2} ) ( \frac{z}{4})^{a-b+ \frac{1}{2}} \sum_{n=0}^{\infty} \frac{(2b-2a-1)_n(b-2a)_n(b-a-\frac{1}{2}+n)}{n!(b)_n} (-1)^n I_{b-a-\frac{1}{2}+n}(\frac{z}{2})$$]
[equation hypergeometric_1f1_3]
=======
[:[/\large $$_1F_1(a, b; z) = e^{ \frac{z}{2} } \Gamma(b-a- \frac{1}{2} ) ( \frac{z}{4})^{a-b+ \frac{1}{2}} \sum_{n=0}^{\infty} \frac{(2b-2a-1)_n(b-2a)_n(b-a-\frac{1}{2}+n)}{n!(b)_n} (-1)^n I_{b-a-\frac{1}{2}+n}(\frac{z}{2})$$]
[$../equations/hypergeometric_1f1_3.svg]]
>>>>>>> 75789acd778461428be66b3bbc3ada25e33a042c
which is particularly useful for ['a [cong] b] and ['z > 0], or /a/ \u226a 1, and ['z < 0].
Then we have Tricomi's approximation (given in simplified form in A&S 13.3.7) [link math_toolkit.hypergeometric.hypergeometric_refs (7)]
<<<<<<< HEAD
[/\large $$_1F_1(a, b; z) = \Gamma(b) e^{ \frac{1}{2} z} \sum_{n=0}^{\infty} 2^{-n}z^n A_n(a,b) E_{b-1+n}(z(\frac{b}{2}-a)) $$]
[equation hypergeometric_1f1_4]
@@ -313,10 +347,24 @@ and
[/\large $$(n+1)A_{n+1} = (n+b-1)A_{n-1} - (b-2a)A_{n-2} \quad;\quad n \geq 2 $$]
[equation hypergeometric_1f1_6]
=======
[:[/\large $$_1F_1(a, b; z) = \Gamma(b) e^{ \frac{1}{2} z} \sum_{n=0}^{\infty} 2^{-n}z^n A_n(a,b) E_{b-1+n}(z(\frac{b}{2}-a)) $$]
[$../equations/hypergeometric_1f1_4.svg]]
with
[:[/\large $$A_0(a,b) = 1, A_1(a,b) = 0, A_2(a,b) = \frac{b}{2}, A_3(a,b)= - \frac{1}{3}(b-2a) $$]
[$../equations/hypergeometric_1f1_5.svg]]
and
[:[/\large $$(n+1)A_{n+1} = (n+b-1)A_{n-1} - (b-2a)A_{n-2} \quad;\quad n \geq 2 $$]
[$../equations/hypergeometric_1f1_6.svg]]
>>>>>>> 75789acd778461428be66b3bbc3ada25e33a042c
With ['E[sub v]] defined as:
[/
[:[/
\begin{equation*}
\begin{split}
E_v(z) & = z^{-\frac{1}{2}v} J_v(2 \sqrt{z})\\
@@ -324,12 +372,17 @@ E_v(-z) & = z^{-\frac{1}{2}v} I_v(2 \sqrt{z})\\
E_v(0) & = \frac{1}{\Gamma(v + 1)}
\end{split}
\end{equation*}]
<<<<<<< HEAD
[equation hypergeometric_1f1_7]
=======
[$../equations/hypergeometric_1f1_7.svg]]
>>>>>>> 75789acd778461428be66b3bbc3ada25e33a042c
This approximation is especially effective when ['a < 0].
For /a/ and /b/ both large and positive and somewhat similar in size then an expansion in terms of gamma functions can be used [link math_toolkit.hypergeometric.hypergeometric_refs (6)]:
<<<<<<< HEAD
[/\large $$_1F_1(a, b; x) = \frac{1}{B(a, b-a)} e^x \sum_{k=0}^{\infty} \frac{(1-a)_k}{k!} \frac{\gamma(b-a+k, x)}{x^{b-a+k}} $$]
[equation hypergeometric_1f1_8]
@@ -337,6 +390,15 @@ For /z/ large we have the asymptotic expansion:
[/\large $$_1F_1(a, b; x) \approx \frac{e^x x^{a-b}}{\Gamma(a)} \sum_{k=0}^{\infty} \frac{(1-a)_k(b-a)_k}{k! x^k} $$]
[equation hypergeometric_1f1_9]
=======
[:[/\large $$_1F_1(a, b; x) = \frac{1}{B(a, b-a)} e^x \sum_{k=0}^{\infty} \frac{(1-a)_k}{k!} \frac{\gamma(b-a+k, x)}{x^{b-a+k}} $$]
[$../equations/hypergeometric_1f1_8.svg]]
For /z/ large we have the asymptotic expansion:
[:[/\large $$_1F_1(a, b; x) \approx \frac{e^x x^{a-b}}{\Gamma(a)} \sum_{k=0}^{\infty} \frac{(1-a)_k(b-a)_k}{k! x^k} $$]
[$../equations/hypergeometric_1f1_9.svg]]
>>>>>>> 75789acd778461428be66b3bbc3ada25e33a042c
which is a special case of the gamma function expansion above.
@@ -350,8 +412,13 @@ That effectively completes the "direct" methods used, the remaining areas are tr
These require extreme care in their use, as they often change the direction of stability, or else are not stable at all.
Sometimes this is a well defined and characterized change in direction: for example with /a,b/ and /z/ all positive recurrence on /a/ via
<<<<<<< HEAD
[/\large $$(b-a) _1F_1(a-1, b; z) + (2a-b+z) _1F_1(a, b; z) -a _1F_1(a+1, b; z) = 0 $$]
[equation hypergeometric_1f1_10]
=======
[:[/\large $$(b-a) _1F_1(a-1, b; z) + (2a-b+z) _1F_1(a, b; z) -a _1F_1(a+1, b; z) = 0 $$]
[$../equations/hypergeometric_1f1_10.svg]]
>>>>>>> 75789acd778461428be66b3bbc3ada25e33a042c
abruptly changes from stable forwards recursion to stable backwards if /2a-b+z/ changes sign.
Thus recurrence on /a/, even when [sub 1]/F/[sub 1](/a/+/N/, /b/, /z/) is strictly increasing, needs careful handling as /a/ \u2192 0.
@@ -370,8 +437,13 @@ the continued fractions associated with the recurrence relations to calculate
and then normalize
either by iterating forwards on /b/ until /b > 0/ and comparison with a reference value, or by using the Wronskian
<<<<<<< HEAD
[/\large $$_1F_1(a, b; z) \frac{d}{dz}z^{1-b}_1F_1(1+a-b, 2-b; z) - z^{1-b}_1F_1(1+a-b, 2-b; z)\frac{d}{dz}{_1F_1}(a, b; z) = (1-b)z^{-b}e^z,$$]
[equation hypergeometric_1f1_11]
=======
[:[/\large $$_1F_1(a, b; z) \frac{d}{dz}z^{1-b}_1F_1(1+a-b, 2-b; z) - z^{1-b}_1F_1(1+a-b, 2-b; z)\frac{d}{dz}{_1F_1}(a, b; z) = (1-b)z^{-b}e^z,$$]
[$../equations/hypergeometric_1f1_11.svg]]
>>>>>>> 75789acd778461428be66b3bbc3ada25e33a042c
which is of course the well known Miller's method [link math_toolkit.hypergeometric.hypergeometric_refs (12)].
@@ -462,8 +534,13 @@ are destroyed by cancellation.
The function `hypergeometric_pFq` returns the result of:
<<<<<<< HEAD
[/\Large $$ _pF_q(\{a_1...a_n\}, \{b_1...b_n\}; z) = \sum_{k=0}^{\infty} \frac{\Pi_{j=1}^n(a_j)_n z^k}{\Pi_{j=1}^n(b_j)_n k!} $$]
[equation hypergeometric_pfq_1]
=======
[:[/\Large $$ _pF_q(\{a_1...a_n\}, \{b_1...b_n\}; z) = \sum_{k=0}^{\infty} \frac{\Pi_{j=1}^n(a_j)_n z^k}{\Pi_{j=1}^n(b_j)_n k!} $$]
[$../equations/hypergeometric_pfq_1.svg]]
>>>>>>> 75789acd778461428be66b3bbc3ada25e33a042c
It is most naturally used via initializer lists as in: