R中的三重积分（如何指定域）

Question

我想在 R 中计算三个变量 f(x,y,z) 的函数的三重积分。我正在使用程序包 cubature 和函数 adaptIntegrate()。被积函数仅在某个域中等于 1（x<y<z，否则为 0），我不知道如何指定。我正在尝试该函数的 2 种不同实现，但其中 none 有效：

#First implementation
fxyz <- function(w) {
x <- w[1]
y <- w[2]
z <- w[3]
x*y*z*(x < y)&(y < z)
}

#Second implementation
fxyz <- function(w) {
x <- w[1]
y <- w[2]
z <- w[3]
if(x<y&y<z)
    out<-1
else
    out<-0
out
}

#Computation of integral
library(cubature)
lower <- rep(0,3)
upper <- rep(1, 3)
adaptIntegrate(f=fxyz, lowerLimit=lower, upperLimit=upper, fDim = 3)

知道如何正确指定域吗？

Answer 1

我不知道 cubature 包，但你可以通过重复应用 base R 的 integrate 函数进行一维积分来做到这一点。

f.xyz <- function(x, y, z) ifelse(x < y & y < z, 1, 0)
f.yz <- Vectorize(function(y, z) integrate(f.xyz, 0, 1, y=y, z=z)$value,
                  vectorize.args="y")
f.z <- Vectorize(function(z) integrate(f.yz, 0, 1, z=z)$value,
                 vectorize.args="z")

integrate(f.z, 0, 1)
# 0.1666632 with absolute error < 9.7e-05

您可能想要使用控制参数来设置数字容差；内部整合的小错误可能会变成外部的大错误。

Answer 2

在您的第一个函数中，return 值是错误的。应该是as.numeric(x<=y)*as.numeric(y<=z)。在您的第二个函数中，您还应该使用 <= 而不是 <，否则 `adapIntegrate 将无法正常工作。您还需要指定最大评估次数。试试这个

library(cubature)
lower <- rep(0,3)
upper <- rep(1,3)

# First implementation (modified)
fxyz <- function(w) {
    x <- w[1]
    y <- w[2]
    z <- w[3]
    as.numeric(x <= y)*as.numeric(y <= z)
}

adaptIntegrate(f=fxyz,lowerLimit=lower,upperLimit=upper,doChecking=TRUE,
          maxEval=2000000,absError=10e-5,tol=1e-5)
#$integral
#[1] 0.1664146
#$error
#[1] 0.0001851699
#$functionEvaluations
#[1] 2000031
#$returnCode
#[1] 0

Answer 3

域 0 <= x <= y <= z <= 1 是 "canonical" 单纯形。要在单纯形上积分，请使用 SimplicialCubature 包。

library(SimplicialCubature)
f <- function(x) 1
S <- CanonicalSimplex(3)

> adaptIntegrateSimplex(function(x) 1, S)
$integral
[1] 0.1666667

$estAbsError
[1] 1.666667e-13

$functionEvaluations
[1] 55

$returnCode
[1] 0

$message
[1] "OK"

请注意，在单纯形上对常量函数 f(x)=1 求积分只是得出单纯形的体积，即 1/6。集成对于此示例没有用。

> SimplexVolume(S)
[1] 0.1666667

R中的三重积分（如何指定域）

Triple integral in R (how to specifying the domain)

r

integral