R中的n维积分，极限为函数答案

【问题标题】：n-dimensional integration in R with limits as functionsR中的n维积分，极限为函数
【发布时间】：2015-12-16 04:43:31
【问题描述】：

我知道诸如cubature 之类的R 包能够对许多积分执行数值积分。但是，在我的情况下，积分的限制不是整数，而是其他变量的函数。

例如，假设我有

f(x, y) = 4*(x-y)*(1-y)，其中 0

我想在 x 的范围内计算 y pracma 包吗？）

如果这还没有实现，我会感到惊讶，但在这种情况下，一个算法的链接会很有帮助。

【问题讨论】：

标签： c r integration numerical

【解决方案1】：

如果您只是添加（或者更确切地说乘以）一个项，该项将函数的值设置为在 (y >= x) 的任何 (x,y) 元组处为零，这实际上非常简单

library(pracma)
fun <- function(x, y) 4*(x-y)*(1-y)*(y < x)
integral2(fun, 0, 2, 0, 1, reltol = 1e-10)
#=============
$Q
[1] 2.833363

$error
[1] 2.394377e-11

积分函数不会接受无限界，但由于y 被限制为（0-1）且x 大于y，那么x 的下界也必须为0，即这就是我以这种方式设置限制的原因。

【讨论】：

条件应该是x >= y 但感谢您的回答
我实际上在我使用*(y < x) 的地方发布了答案。然后我看了看它并决定它被翻转并改变了它，但现在我审查了逻辑我同意你是对的并且会恢复。

【解决方案2】：

您可以将域分成三个三角形（制作图片）。通过在列中给出它们的顶点来定义这三个三角形：

T1 <- cbind(c(0,0), c(1,0), c(1,1))
T2 <- cbind(c(1,0), c(2,0), c(1,1))
T3 <- cbind(c(1,1), c(2,1), c(2,0))

然后您可以使用SimplicialCubature 包。先定义一个数组中三个三角形的并集：

Domain <- abind::abind(T1, T2, T3, along=3)

定义被积函数：

f <- function(xy) 4*(xy[1]-xy[2])*(1-xy[2])

然后：

library(SimplicialCubature)
> adaptIntegrateSimplex(f, Domain)
$integral
[1] 2.833333

$estAbsError
[1] 2.833333e-12

$functionEvaluations
[1] 96

积分的准确值为17/6 = 2.833333....。因此，我们得到的近似值优于 pracma 在另一个答案中给出的近似值（不足为奇：带有术语 (y<x) 的被积函数不平滑）。

对于这个例子，我们甚至可以使用SimplicialCubature 得到更好的结果。被积函数是多项式：f(x,y) = 4*x - 4*xy - 4*y + 4*y²。 SimplicialCubature 包有一个精确的多项式方法。

> p <- definePoly(c(4,-4,-4,4), rbind(c(1,0),c(1,1),c(0,1),c(0,2)))
> printPoly(p)
Polynomial has 4 terms and 2 variables
4 * x[1]^1       (degree= 1 )
  -4 * x[1]^1 * x[2]^1       (degree= 2 )
  -4 * x[2]^1       (degree= 1 )
  + 4 * x[2]^2       (degree= 2 )
> integrateSimplexPolynomial(p, Domain)
$integral
[1] 2.833333

$functionEvaluations
[1] 9

【讨论】：