puis-je faire confiance à cette triple intégrale numérique de Matlab?

Gens des sciences informatiques:

J'ai initialement posté cette question sur Math Stack Exchange et quelqu'un a commenté que je pourrais obtenir des réponses "bien meilleures" ici:

Je suis novice en méthodes numériques et en Matlab. J'essaie d'évaluer la somme suivante de deux triples intégrales (elle peut évidemment être écrite plus simplement, mais vous ne pouvez toujours pas l'évaluer symboliquement (?)). J'ai du mal à faire fonctionner le ici, donc je l'ai divisé à contrecœur en morceaux ici: je veux trouver la somme de $\LaTeX$

\frac{2}{((1 / 0.3) - 1)^{2}} (\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}) d t d r_{0} d r_{1}),

$\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),$

\frac{2}{((1 / 0.3) - 1)^{2}} (\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}) d t d r_{0} d r_{1}),

$\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),$

où

F_{1} (r_{0}, r_{1}, t) = \frac{t^{2} r_{0}^{3} * (0.3)^{3}}{2 r_{1}^{3} \sqrt{π}}

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

F_{2} (r_{0}, r_{1}, t) = \frac{(0.3)^{3} π^{3 / 2} (r_{0} + r_{1} - t)^{4} (t^{2} + 2 t (r_{0} + r_{1}) - 3 (r_{1} - r_{0})^{2})^{2}}{288 (\frac{4}{3} π r_{0}^{3}) (\frac{4}{3} π r_{1}^{3})} .

$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 mars 2013): Quelqu'un a répondu qu'il avait demandé à Mathematica de faire symboliquement les intégrales. J'ai juste essayé de le faire (avec des versions simplifiées des intégrales) et Mathematica ne pouvait que faire les deux externes du premier, et a calé sur le second. J'apprécierais de l'aide. Voici ce que j'ai fait.:

J'ai essayé d'évaluer

\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}} ré t ré r_{1} ré r_{2}

$\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

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

et Mathematica revient (j'ai eu des problèmes avec le ici parce que le résultat est long. Je l'ai divisé en deux équations. Si quelqu'un connaît un bon moyen d'afficher cela, dites-le moi): $\LaTeX$

\int_{1}^{2} \frac{1}{64 r 2^{2}} e^{- 1 - r 2^{2}} (2 e^{2 r 2} (25 + r 2 (19 + 2 r 2 (1 + r 2))) -

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

e^{1 + r 2^{2}} (32 r 2 (2 + r 2^{2})) + \sqrt{π} (11 + 4 r 2^{2} (9 + r 2^{2})) Erf [1 - r 2]) ré r 2.

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

Ensuite, j'ai essayé d'évaluer

\int_{1}^{2} \int_{1}^{r_{2}} \int_{r_{2} - r_{1}}^{r_{2} + r_{1}} \dots

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

\dots \frac{\exp (- t^{2}) (r_{1} + r_{2} - t)^{4} (t^{2} + 2 t (r_{1} + r_{2}) - 3 (r_{2} - r_{1})^{2})^{2}}{r_{1}^{3} r_{2}^{3}} ré t ré r_{} ré r_{2}

$\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$

en utilisant

Intégrer [(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}]

juste maintenant, et Mathematica n'a pas retourné de réponse après environ une demi-heure (mais j'ai des problèmes de réseau informatique en ce moment, et ils peuvent être à blâmer).

[FIN DE LA MODIFICATION DU 2 MARS]

J'ai utilisé la commande "triplequad" de Matlab, sans options supplémentaires. J'ai géré les limites variables de l'intégration au moyen de fonctions heaviside, car je ne connaissais pas d'autre moyen de le faire. Matlab m'a donné . $0.007164820144202$

Je sais que Matlab est un bon logiciel, mais j'ai entendu dire que les intégrales triples numériques sont difficiles à faire avec précision et que les mathématiciens sont censés être sceptiques, donc je veux un moyen de vérifier l'exactitude de cette réponse. Les intégrales donnent la valeur attendue d'une certaine expérience (si quelqu'un le souhaite, je peux éditer cette question pour décrire l'expérience): J'ai implémenté l'expérience dans Matlab en utilisant des nombres générés de manière appropriée et aléatoire, un million de fois, et en moyenne les résultats. J'ai répété ce processus quatre fois. Voici les résultats (je m'excuse si j'ai mal utilisé le mot «procès»):

Essai 1: $0.007133292603256$

Essai 2: $0.007120455071989$

Essai 3: $0.007062595022049$

Essai 4: $0.007154940168452$

Essai 5: $0.007215000289130$

Bien que chaque essai ait utilisé un million d'échantillons, les valeurs de simulation ne concordent que dans le premier chiffre significatif. Ils ne sont pas assez proches les uns des autres pour que je puisse déterminer si la triple intégrale numérique est exacte.

Alors, quelqu'un peut-il me dire si je peux faire confiance au résultat de "triplequad" ici, et dans quelles circonstances peut-on lui faire confiance en général?

Une suggestion que j'ai eue à Math Stack Exchange était d'essayer d'autres logiciels comme Mathematica, Octave, Maple et SciPy. Est-ce un bon conseil? Est-ce que les gens font du travail numérique dans Mathematica et Maple? Octave est une sorte de clone Matlab, alors puis-je supposer qu'il utilise les mêmes algorithmes d'intégration? Je n'ai même jamais entendu parler de SciPy auparavant et j'aimerais avoir des opinions à ce sujet.

MISE À JOUR: Quelqu'un de Math Stack Exchange l'a fait dans Maple et a obtenu . Cela correspond à trois chiffres importants. C'est un bon signe. $0.007163085468$

En outre, j'apprécierais des suggestions sur la façon d'entrer une expression longue et multiligne dans dans Stack Exchange. Pouvez-vous utiliser l'environnement "aligné" ici? J'ai essayé et je n'ai pas réussi à le faire fonctionner. $\LaTeX$

matlab Stefan Smith
la source

Vos résultats de simulation sont parfaitement cohérents avec la valeur numérique renvoyée par Matlab: leur moyenne de n'est que de erreur standard inférieure à ce que Matlab a renvoyé. FWIW, Mathematica renvoie . Il peut également évaluer symboliquement ces intégrales en termes de polynômes et de fonctions d'erreur.

0.00713726

$0.00713726$

- 1.11

$-1.11$

0.00716308537 \dots

$0.00716308537\ldots$

whuber

@whuber Merci. Je pourrais jurer que je l'ai essayé symboliquement dans Maple et que Maple n'a pas pu le faire. Je vais réessayer dans Maple, et si ça ne marche pas, je vais essayer dans Mathematica. BTW, j'ai fait une intégrale similaire dans Maple, et j'ai obtenu une énorme réponse symbolique. Cela semblait être une somme et une différence de très grands nombres dont le grand total était assez petit. Je soupçonne qu'une erreur d'arrondi était probable dans la réponse finale. Dans un problème comme celui-ci, devez-vous utiliser la réponse symbolique, ou simplement faire l'intégrale numériquement?

Stefan Smith

Les réponses symboliques ont l'avantage d'être des combinaisons de fonctions qui (souvent) peuvent être calculées efficacement avec une précision arbitraire. Habituellement aussi, la solution symbolique se prête également à un recalcul rapide lorsque les paramètres varient. Pour ces raisons, il vaut souvent la peine de rechercher une solution symbolique.

whuber

@whuber: J'ai essayé de faire des intégrales essentiellement équivalentes (en modifiant certaines des constantes et en supprimant certaines constantes multiplicatives) dans Mathematica, et Mathematica ne pouvait faire que les deux intégrations externes de la première intégrale, et semble avoir calé sur la seconde. J'ai posté mon code et mes résultats ci-dessus.

Stefan Smith

Concernant l'édition du 2 mars: En réduisant symboliquement l'intégrale triple à une seule intégrale (dans la première moitié de vos intégrales), vous avez accompli beaucoup de choses. L'intégrande se comporte très bien et peut être intégrée numériquement avec une précision extrêmement élevée en une fraction de seconde.

whuber

puis-je faire confiance à cette triple intégrale numérique de Matlab?

Réponses: