R中的双积分实现答案

【问题标题】：Double Integral implementation in RR中的双积分实现
【发布时间】：2019-12-05 00:12:39
【问题描述】：

我正在使用 R 中的 mosaicCalc 包处理双积分。我无法获得第二个双积分的正确结果。

这是产生正确结果 (pi) 的第一个双积分的代码。

one = makeFun(1 ~ y + x)
bx.y = antiD(one(y = y, x = x) ~ y)
bx.yx = antiD(bx.y(y = 1+cos(x), x = 2) ~ x)
bx.yx(x = pi) - bx.yx(x = 0)
# [1] 3.141593

这是第二个双积分，根据 Wolfram 的说法，它的正确结果应该是 0.684853

one = makeFun(1 ~ y + x)
bx.y = antiD(one(y = y, x = x) ~ y)
bx.yx = antiD(bx.y(y = 1/2, x = sin(x)) ~ x)
bx.yx(x = 5*pi/6) - bx.yx(x = pi/6)
# [2] 1.047198

【问题讨论】：

欢迎来到 SO，Dessel！这个问题缺乏很多信息让我们提供任何真正的帮助。请让这个问题可重现。这包括示例代码（包括列出非基础 R 包）、示例明确数据（例如，dput(head(x)) 或data.frame(x=...,y=...)）和预期输出。参考：stackoverflow.com/questions/5963269、stackoverflow.com/help/mcve 和 stackoverflow.com/tags/r/info。
现在您已经添加了图片......到目前为止您尝试了什么？（说“没有成功”对我们没有多大帮助。如果您遇到错误，请提供代码和错误，我们可以提供帮助。如果您收到警告和/或错误结果，我们需要代码，实际输出，以及为什么不正确。）
我用正确的积分和代码更新了我的问题。对于给您带来的不便，我深表歉意。
好多了，谢谢 Dessel！

标签： r double integral

【解决方案1】：

首先我需要说服自己 Mathematica 是正确的。是的，我想这源于我对权威的不信任，但这并不难。需要认识到从上到下的统一积分就是上减下的差，所以可以简化为单变量问题，用R的integrate函数解决：

integrate( function(x){sin(x)-0.5},lower=pi/6,upper=5*pi/6)
0.6848533 with absolute error < 7.6e-15

因此，我在“mosaicCalc”框架中尝试了相同的策略：

> one = makeFun(1 ~ y + x)
> bx.y = antiD(one(y = y, x = x) ~ y)
> bx.yx = antiD(bx.y(y = sin(x)-1/2, x = x) ~ x)
> bx.yx(x = 5*pi/6) - bx.yx(x = pi/6)
[1] 0.6848533

这得到了正确的答案，但它似乎没有正确地表示和保留“功能级联”（发明了一个我以前从未听说过的术语）。我希望限制以反映更通用的调用功能集的方式出现，所以最终想出了这个看起来令人满意的方法：

one = makeFun(1 ~ y + x)
bx.y = antiD(one(y = y, x = x) ~ y)
bx.yx = antiD(bx.y(x=x, y = sin(x)) - bx.y(x=x,y=1/2) ~ x)
bx.yx(x = 5*pi/6) - bx.yx(x = pi/6)
[1] 0.6848533

这确实支持更复杂的函数被积函数。我没有设置 Mathematica 来与任何东西进行集成，但使用上面的镶嵌计算器设置，我得到 1.284286 for x^2 作为被积函数。你可能想检查一下。

在解决这个问题时，在我看来，整合的顺序是颠倒的。在我从 40 年前（不是 50 年前）对这些问题的模糊记忆中，似乎我一直使用 dx 计算作为内部计算，但我确实意识到这是任意的。无论如何，您在第二个反导数中分配给 x 和 y 值的角色似乎并不合适。您得到了从 x=pi/6 和 y=0.5 的下限到 x=5*pi/6, y=1) 的上限的 2D 积分结果

library(cubature) # package capable of 2D-integration with fixed limits
adaptIntegrate(function(x){1}, lower=c(a=pi/6, b=0.5), upper=c(a=5*pi/6,b=1 ))
$integral
[1] 1.047198

$error
[1] 0

$functionEvaluations
[1] 17

$returnCode
[1] 0

【讨论】：