复杂函数与 Cubature C 包的数值积分答案

【问题标题】：Numerical integration of a complex function with Cubature C package复杂函数与 Cubature C 包的数值积分
【发布时间】：2013-05-28 14:22:36
【问题描述】：

我想使用Cubature C 包来执行复杂函数的多维积分。对于正方形 [0,1]x[0,1] 上的非常简单的函数 f(x,y) = x + y*i，我尝试按以下方式执行此操作。确切的结果是 0.5 + 0.5i。

#include <stdio.h>
#include <math.h>
#include <complex.h>
#include "../cubature.h"

int f(unsigned ndim, const double *x, void *fdata, unsigned fdim, double *fval);

int main(void)
{
    double xmin[2] = {0,0}, xmax[2] = {1,1}, val, err;
    hcubature(1, f, NULL, 2, xmin, xmax, 0, 0, 1e-3, ERROR_PAIRED, &val, &err);
    printf("Computed integral = %f+%fi +/- %f\n", creal(val),cimag(val), err);
}

int f(unsigned ndim, const double *x, void *fdata, unsigned fdim, double *fval) 
{
    fval[0] = x[0]+x[1]*I;
    return 0;
}

但我得到的是Computed integral = 0.500000+0.000000i +/- 0.000000，即没有出现虚部。如果我放一个纯虚数被积函数（例如 x*i），我总是得到 0 作为结果。

我做错了什么？你知道用 C 计算复杂函数积分的更好方法吗？

【问题讨论】：

我对这个包一无所知，但我不清楚fval[0] = x[0]+x[1]*I;如何成功地将一个复数打包成一个单独的double
R 标签是真的吗？ R（语言）中有一个 cubature 包。
是的，有！你可以找到它here
@AakashM 你是绝对正确的。
@RiccardoCavallari 我知道，但我的意思是：这个问题是否与 R 包有关？

标签： c math numerical-methods complex-numbers

【解决方案1】：

2 个问题：

1) 你声明double val，但你希望它是double val[2]。您在 cimag(val) 中作为 double 的原始 val 将始终返回 0.0。改为

double xmin[2] = {0,0}, xmax[2] = {1,1}, val[2] /* not val */, err;
hcubature(2 /* not 1 */, f, NULL, 2, xmin, xmax, 0, 0, 1e-3, ERROR_PAIRED, /* & */val, &err);
printf("Computed integral = %f+%fi +/- %f\n", val[0], val[1], err);

2) 在f() 中，您似乎没有正确计算f(x,y) = x + y*i。 fval[0] 是双精度数，无法保存您尝试分配的复杂答案。

int f(unsigned ndim, const double *x, void *fdata, unsigned fdim, double *fval) {
    fval[0] = x[0];
    fval[1] = x[1];
    return 0;
}

抱歉，我现在没有要测试的“../cubature.h”。

【讨论】：

【解决方案2】：

定义 val 时，请执行以下操作：

double *val=(double *) malloc(sizeof(double) * integrand_fdim);

其中integrand_fdim 等于unsigned fdim - 你的被积函数f 的第一个参数。

P.S.：@chux 的回答产生了 Bus 10 或 segfault...

【讨论】：