将Mathematica积分计算转换为R语言时的区域限制问题求助

原Mathematica代码逻辑 #

原代码是在圆心(0,0)、半径10的圆盘区域内，对满足|x - xprime| < dimx/2且|y - yprime| < dimy/2的矩形子区域，积分函数1/(20π√(x²+y²))。

现有R代码的问题 #

• 圆盘半径错误：原代码为10，误写为5
• 函数系数错误：原函数分母是20π，误写为10π
• 积分上下界顺序错误：calculus::integral要求下限在前、上限在后，现有代码顺序颠倒
• 缺少区域交集判断：未同时验证点是否在矩形和圆盘的交集内，不符合原Boole条件的逻辑

修正后的R代码 #

library(calculus)

# 与原Mathematica一致的参数
xprime <- 0.5
yprime <- 3
dimx <- 20/1000
dimy <- 20/1000
disk_radius <- 10

# 定义积分函数：仅在矩形与圆盘的交集区域计算函数值
pb_function <- function(x, y) {
  in_rect <- abs(x - xprime) < dimx/2 & abs(y - yprime) < dimy/2
  in_disk <- sqrt(x^2 + y^2) <= disk_radius
  ifelse(in_rect & in_disk, 1/(20 * pi * sqrt(x^2 + y^2)), 0)
}

# 执行积分：修正上下界顺序
result <- integral(
  pb_function,
  list(
    x = c(xprime - dimx/2, xprime + dimx/2),
    y = c(yprime - dimy/2, yprime + dimy/2)
  )
)

print(result)

说明 #

因为目标矩形的尺寸很小（仅20mm），远小于圆盘半径10，所以直接以矩形为积分范围，在函数内过滤圆盘外的点，既符合原逻辑又保证计算效率。

内容的提问来源于stack exchange，提问作者J. W. Rozelle