Timeline for Numerical integration using interval arithmetic, nowadays
Current License: CC BY-SA 3.0
14 events
| when toggle format | what | by | license | comment | |
|---|---|---|---|---|---|
| Sep 8, 2016 at 5:09 | history | bounty ended | CommunityBot | ||
| Sep 7, 2016 at 15:08 | vote | accept | H A Helfgott | ||
| Sep 5, 2016 at 22:40 | comment | added | Ingo Blechschmidt | Dear Fredrik, let me add that I think that your work on Arb is extremely valuable! A free (as in freedom) and curated software library for interval arithmetic is very useful. Also I've been enjoying your insightful blog posts. | |
| Sep 1, 2016 at 7:04 | comment | added | Fredrik Johansson | The build steps used in .build_dependencies and .travis.yml might be helpful. Send me an email if you're still having problems with it. | |
| Sep 1, 2016 at 3:09 | comment | added | H A Helfgott | Thank you very much! Of course, first I must succeed in getting anything to compile - somehow the routines in Flint refuse to be found at the linking stage... | |
| Sep 1, 2016 at 2:38 | comment | added | Fredrik Johansson | Certainly. Here is a complete .c file implementing the composite trapezoidal rule. Let me know if you run into any issues .gist.github.com/fredrik-johansson/… | |
| Sep 1, 2016 at 0:57 | comment | added | H A Helfgott | If it is easy to implement, might it be possible to ask for a tutorial on how to do it in arb? | |
| Sep 1, 2016 at 0:02 | comment | added | Fredrik Johansson | The best way to use the derivatives also depends heavily on the function. You might get tight enclosures for $f$ but not for $f^{(d)}$, or vice versa. Evaluation of $f$ might be fast while $f'$, $f''$, ... are slow, or $f$ might be relatively expensive while the first few derivatives are almost free in comparison (this is the case with $f = \operatorname{erf}$). | |
| Aug 31, 2016 at 23:43 | comment | added | Fredrik Johansson |
In Arb, I'm fairly systematically implementing special functions with support for truncated power series input and output (= automatic differentiation), so you can get arbitrary derivatives easily. Again, this isn't yet exposed in Sage, though it might be of some use to note that I have a separate Python interface that does wrap the power series conveniently. For example, you can do ctx.cap=3; x=arb_series([0,1]); (1+x*x).sqrt().erf() to produce $a_0 + a_1 x + a_2 x^2 + O(x^3) = \operatorname{erf}(\sqrt{1 + x^2})$.
|
|
| Aug 31, 2016 at 23:33 | comment | added | Fredrik Johansson | You're right, and that's indeed easy to implement. The drawback is that for $(d+1)$th order error bounds, you need to bound $|f^{(d)}|$ on intervals, which can be difficult when $d$ is large. With the complex method, you only need to bound $|f|$ itself on intervals; the bounds can then be reused for any $d$. (Taylor methods avoid evaluation on intervals completely.) All of this only really matters if you need so high precision that say $d = 2$ or $d = 4$ would require too large $N$; for modest precision, it would make sense most of the time to use a low order method and increase $N$. | |
| Aug 31, 2016 at 20:05 | comment | added | H A Helfgott | Of course I can do this by differentiating $f$ by hand (or in SAGE) and then typing in the resulting expression (very hairy, in what I am doing) - however, this is error-prone, and the procedure ought to be automated. think?What do you | |
| Aug 31, 2016 at 20:04 | comment | added | H A Helfgott | Let us work with the simplest kind of numerical integration - higher-order versions should work in the same way. We want to approximate the integral by a sum of type $\sum_{n=1}^N f(a_i)$, where $a_i = a + i (b-a)/N$. The error term is $\leq ((b-a)/N)^2 \sum_{n=1}^N \max_{a_{i-1}\leq x\leq a_i} |f'(x)|$. Now, we can compute the maximum within the sum simply by using interval arithmetic to bound the image of [a_{i-1},a_i] under $f'$. All we need, then, is interval arithmetic and symbolic differentiation. | |
| Aug 31, 2016 at 19:51 | comment | added | H A Helfgott | Thanks - I'll try this. Still, I do not quite get why the following simpler strategy (which requires no use of complex analysis) couldn't be automated. | |
| Aug 31, 2016 at 8:35 | history | answered | Fredrik Johansson | CC BY-SA 3.0 |