dapatkah saya mempercayai integral tiga angka ini dari Matlab?

Orang-orang Ilmu Komputasi:

Saya awalnya memposting pertanyaan ini di Math Stack Exchange dan seseorang berkomentar bahwa saya mungkin mendapatkan jawaban "jauh lebih baik" di sini:

Saya seorang pemula dalam metode numerik dan Matlab. Saya mencoba untuk mengevaluasi jumlah berikut dari dua integral tiga (ini jelas dapat ditulis lebih sederhana, tetapi Anda masih tidak dapat mengevaluasinya secara simbolis (?)). Saya kesulitan mendapatkan bekerja di sini, jadi saya dengan enggan memecah-mecahnya menjadi beberapa bagian di sini: Saya ingin mencari jumlah$\LaTeX$

$$\frac{2}{((1/0.3) - 1)^2}\left(\int_1^{1/0.3}\int_1^{r_1}\int_0^{r_1-r_0}F_1(r_0,r_1,t)\exp(-\frac{(0.3)^2 t^2}{4})\,dt\,dr_0\,dr_1 \right),$$

dan

$$\frac{2}{((1/0.3) - 1)^2}\left(\int_1^{1/0.3}\int_1^{r_1}\int_{r_1-r_0}^{r_1+r_0} F_2(r_0,r_1,t)\exp(-\frac{(0.3)^2 t^2}{4})\,dt\,dr_0\,dr_1 \right),$$

dimana

$$F_1(r_0,r_1,t)=\frac{t^2 r_0^3*(0.3)^3}{2r_1^3\sqrt{\pi}}$$

and

$$F_2(r_0,r_1,t)=\frac{(0.3)^3\pi^{3/2}(r_0+r_1-t)^4 (t^2+2t(r_0+r_1)-3(r_1-r_0)^2)^2}{288(\frac{4}{3}\pi r_0^3)(\frac{4}{3}\pi r_1^3)}.$$

EDIT (2 Maret 2013): Seseorang menjawab bahwa mereka mendapat Mathematica untuk melakukan integral secara simbolis. Saya hanya berusaha melakukan ini (dengan versi yang disederhanakan dari integral) dan Mathematica hanya bisa melakukan dua yang terluar dari yang pertama, dan berhenti pada yang kedua. Saya akan sangat menghargai bantuan. Inilah yang saya lakukan .:

Saya berusaha mengevaluasi

$$\int_1^2 \int_1^{r_2} \int_0^{r_2-r_1} \frac{r_1^3 t^2 \exp(-t^2)}{r_2^3}\,dt\,dr_1\,dr_2$$ via

Integrate[r1^3/r2^3*t^2*Exp (-t^2), {t, 0, r2 - r1}, {r1, 1, r2}, {r2, 1, 2}]

and Mathematica returns (I had trouble with the $\LaTeX$ here because the result is long. I broke it into two equations. if anyone knows a good way to display this please tell me):

$$\int_1^2 \frac{1}{64r2^2} e^{-1-r2^2}(2e^{2r2}(25+r2(19+2r2(1+r2)))-$$

$$e^{1+r2^2}(32r2(2+r2^2)) +\sqrt{\pi}(11+4r2^2(9+r2^2))\operatorname{Erf}[1-r2])\,dr2.$$

Then I tried to evaluate

$$\int_1^2\int_1^{r_2}\int_{r_2-r_1}^{r_2+r_1} \ldots \qquad \qquad \qquad$$

$$\ldots\frac{\exp(-t^2)(r_1+r_2-t)^4(t^2+2t(r_1+r_2)-3(r_2-r_1)^2)^2}{r_1^3 r_2^3}\,dt\,dr_\,dr_2$$

using

Integrate[(r1 + r2 - t)^4*(t^2 + 2*t*(r1 + r2) - 3*(r2 - r1)^2)^2* Exp[-t^2]/r1^3/r2^3, {r2, 1, 2}, {r1, 1, r2}, {t, r2-r1, r2 + r1}]

just now, and Mathematica has not returned an answer after about half an hour (but I am having computer network problems right now, and they may be to blame).

[END OF MARCH 2 EDIT]

I used Matlab's "triplequad" command, with no extra options. I handled the variable limits of integration by means of heaviside functions, because I didn't know any other way to do it. Matlab gave me $0.007164820144202$.

I know Matlab is good software, but I have heard that numerical triple integrals are hard to do accurately, and mathematicians are supposed to be skeptical, so I want some way to verify the accuracy of this answer. The integrals give the expected value of a certain experiment (if anyone wants, I can edit this question to describe the experiment): I implemented the experiment in Matlab using appropriately randomly generated numbers, a million times, and averaged the results. I repeated this process four times. Here are the results (I apologize if I have used the word "trial" improperly):

Trial 1: $0.007133292603256$

Trial 2: $0.007120455071989$

Trial 3: $0.007062595022049$

Trial 4: $0.007154940168452$

Trial 5: $0.007215000289130$

Although each trial used a million samples, the simulation values only agree in the first significant digit. They are not close enough to each to each other for me to determine whether the numerical triple integral is accurate.

So can anyone tell me whether I can trust the result of "triplequad" here, and under what circumstances one can trust it in general?

One suggestion I got at Math Stack Exchange was to try other software like Mathematica, Octave, Maple, and SciPy. Is this good advice? Do people actually do numerical work in Mathematica and Maple? Octave is kind of a Matlab clone, so can I assume it uses the same integration algorithms? I haven't even heard of SciPy before and would appreciate any opinions about it.

---

UPDATE: Someone from Math Stack Exchange did it in Maple and got $0.007163085468$. That is agreement to three significant figures. That is a good sign.

Also, I would appreciate suggestions on how to enter long, multi-line expression in $\LaTeX$ in Stack Exchange. Can you use the "aligned" environment here? I tried, and I couldn't get it to work.

First of all, it is not the software (or at least it shouldn't be) that determines the quality of the solution to a problem, it's the quality and appropriateness of the algorithm that is applied. You should check what algorithm is being used by triplequad in Matlab (I would guess it uses a nested adaptive Gaussian quadrature). And you should check what the requested tolerances are (required absolute and relative tolerance). Chances are that, by default, it only asks for $10^{-8}$ relative precision.

The answer coming from Maple is probably done by Computer Algebra and maybe it could find a closed solution which was then evaluated using double precision floating point. This has the advantage that you are not approximating the integral by a finite summation (and hence are introducing approximation errors) but the Computer Algebra System will find an expression for the integral which can then be evaluated. Of course, care must be taken when evaluating this expression (for round-off).

If you would want to do this with SciPy, you would also need to resort to nested adaptive Gaussian quadrature using the underlying Quadpack (Piessens et al.) routines. In Octave, you'll have the same approach. And I wouldn't be too surprised if Matlab also uses Quadpack as quadrature engine (since it is the reference).