将数组传递给函数指针(在 c 中)?

Passing array into function pointer (in c)?

我想做的是使用 Runge Kutta 求解二阶 ODE。我有 RK 方法的轮廓,并且 ODE 本身已建立,但我无法继续进行。我有一个数组中的 ODE-s,我尝试将这个数组传递给 RK,但是在 f(dy[i]) 部分它给出了 "expected double* but type double instead" 错误。当我尝试将两个 ODE-s 分成两个函数,并将它们传递给 RK 时,它只适用于第一个函数。传递 dy 的元素而不是函数指针不是一种选择,因为 RK 本身应该能够求解包含任意数量变量的 ODE。

如何正确添加数组并使用它? 这是我对应的一段代码:

double* harmOsc(double* p, double t, double* y, double* dy, int n)
{
    double k = p[0];
    double m = p[1];
    double x = y[0];
    double v = y[1];
    dy[0] = v;
    dy[1] = -k/m*x;
    return(dy);
}


double* HarmOsc1(double* dy, double t)
{
    return(dy);
}


void RK4(
    double t,       //independent variable
    double dt,      //stepsize

    double *y,      //variables
    double *dy,     //derivatives
    int n,          //number of equations
    double* (*f)(double*, double))
{

    int j;
    double* l = (double*)calloc(4*n,sizeof(double));
    for(j=0;j<n;j++)
    {
        l[0*n+j] = dt*f(dy[j],t);
        l[1*n+j] = dt*f(dy[j]+0.5*l[0*n+j],t+0.5*dt);
        l[2*n+j] = dt*f(dy[j]+0.5*l[1*n+j],t+0.5*dt);
        l[3*n+j] = dt*f(dy[j]+0.5*l[2*n+j],t+0.5*dt);
        y[j] = y[j] + (l[0*n+j] + 2*l[1*n+j] + 2*l[2*n+j] + l[3*n+j])/6;
    }
}

int main(int argc, char *argv[])
{
    double* p = (double*)calloc(2,sizeof(double));
    p[0] = 15; p[1] = 140;
    double* y = (double*)calloc(2,sizeof(double));
    y[0] = 12.4; y[1] = 1.1;
    double t=0;

    double* dy = (double*)calloc(2,sizeof(double));
    dy = harmOsc(p,t,y,dy,2);
    dy = HarmOsc1(dy,t);

    RK4(t,0.05,p,y,dy,2,&HarmOsc1);
    printf("%f, %f\n",dy[1],dy[0]);
}

HarmOsc1是我调用的函数,所以它有需要的参数量。

当然还有我收到的警告:

RK3.c:55:13: error: incompatible type for argument 1 of ‘f’

RK3.c:55:13: note: expected ‘double *’ but argument is of type ‘double’

RK3.c:56:6: error: incompatible type for argument 1 of ‘f’

RK3.c:56:6: note: expected ‘double *’ but argument is of type ‘double’

RK3.c:57:6: error: incompatible type for argument 1 of ‘f’

RK3.c:57:6: note: expected ‘double *’ but argument is of type ‘double’

RK3.c:58:6: error: incompatible type for argument 1 of ‘f’

RK3.c:58:6: note: expected ‘double *’ but argument is of type ‘double’

函数 f 需要一个 double *,即它的第一个参数指向 double 的指针。但是,在调用 f 的每个地方,您传递的是 double 值,而不是指向 double 的指针,也不是数组(衰减为指向第一个值的指针):

    // here ------------v
    l[0*n+j] = dt*f(  dy[j],                 t);
    l[1*n+j] = dt*f(  dy[j]+0.5*l[0*n+j],    t+0.5*dt);
    l[2*n+j] = dt*f(  dy[j]+0.5*l[1*n+j],    t+0.5*dt);

如果此函数的第一个参数需要是一个数组,则向其传递一个数组,而不是单个值。

f 期望第一个参数是 double*,而你用 dy[j] 调用 f,它是一个 double。 只需将 dy[j] 更改为 dy + j 即可。

您可能需要一个临时变量

double tempv;

然后在循环里面,

tempv = dy[j];
...
tempv = dy[j] + ...;

用tempv的地址调用f,

f(&tempv,...)

根据评论中提供的说明,我们确定您的函数 RK4() 应该使用 Runge-Kutta 方法对形式为

的多个 ODE 中的每一个的一个值进行数值求解
y'(t) = f(y(t), t)

具有以下形式的初始值条件

y(t0) = y0

其中 t 是标量,每个 y 都是标量值,并且 f() 对于所有 ODE 都是相同的(或者至少由相同的 C 表示功能)。

我们进一步确定函数参数具有以下含义:

t    equivalent to t0 above
dt   the distance from t to the point at which the solutions are to be computed
y    points to an array in which to return the solutions
dy   points to the initial function values (y0 above) for all the ODEs
n    the number of ODEs to solve
f    the function `f()` above

你有各种各样的问题,但我认为大多数问题最终都是从你对 RK4() 的参数 f 的声明开始的。请特别注意,ODE 的上述形式要求 f() 接受一个与每个 y 的值具有相同维数的值和另一个与每个 y 具有相同维数/数字的值作为参数] 的参数,以及 return 与每个 y 的值具有相同维度的值。但是我们已经确定 y 是每个 标量 函数 一个标量 参数,所以 f 应该接受两个doubles 和 return 一个 double。这将使我们得到这个声明:

void RK4(
    double t,       //initial point
    double dt,      //delta
    double *y,      //result values
    double *dy,     //initial values
    int n,          //number of equations
    double (*f)(double, double)) {

您不需要 f 在一个 运行 中计算多个 ODE 的结果,因为您的 RK4() 遍历 ODE 并分别为每个 ODE 调用 f() .

接下来,让我们看看您的变量 l。您正在为每个输入 ODE 动态分配足够的 space 来存储所有四个 RK 常量(而不是释放它),但这是无用的。每组 RK 常数仅使用一次,因此在继续下一个方程式时无需记住之前的常数。因此,您只需要 space 四个常量,因为这是一个固定数字,您不需要动态分配。您甚至不会真正从使用数组中获得任何好处;我只使用四个标量变量,在循环体内声明。

    int j;
    // double* l = (double*)calloc(4*n,sizeof(double));
    for (j = 0; j < n; j++) {
        double k1 = dt * f(dy[j], t);
        double k2 = dt * f(dy[j] + 0.5 * k1, t + 0.5 * dt);
        double k3 = dt * f(dy[j] + 0.5 * k2, t + 0.5 * dt);
        double k4 = dt * f(dy[j] + k3, t + dt);  // <-- note corrections

观察现在声明 f 的参数和 return 类型符合 RK 方程对输入 ODE 形式的要求。

最后,计算结果:

        y[j] = dy[j] + (k1 + 2 * k2 + 2 * k3 + k4) / 6;

请注意最后一行与您的版本之间的本质区别 -- RK 计算增量 y,但您的初始 y 值在 dy 中,而不是y.

就是这样:

    }
}

此时我观察到,我在解决您的问题时遇到的最大问题是不了解您尝试做的事情的细节。这是由几件事引起的,其中很大

  • 您的变量和参数名称简短且缺乏表达力,更糟糕的是,它们 似乎 有某种意义,其中有几个实际上代表了不同的东西比他们的名字所暗示的要多

  • 代码文档很少,而且显然不正确。这似乎与前一点有点相似,因为根据有限(但不完全正确)的参数文档,函数参数命名更有意义。

  • 代码本身的帮助有限,因为问题本身存在缺陷。

带回家:编写好的代码文档。使用完整的句子。描述每个函数对每个参数的期望,以及它承诺如何处理这些参数。描述 return 值。记录函数执行的任何错误处理。当您执行此操作时,尝试像想要使用您的功能但无法查阅其实现以了解其功能的人一样思考——此人需要了解什么?事实上,我建议在编写每个函数的实现之前为每个函数编写文档。当您发现需要时,请务必更新文档,但文档可以帮助您保持正轨,并且先编写它们可以帮助您写作,就好像您看不到实现一样。

在 C 中,不建议使用 return 数组。这是可能的,但管理负责释放此内存的代码很麻烦。此外,它违背了一些生成高性能代码的 C 哲学,不必要地重复生成本质上相同的数组,只是在不久之后将其销毁。您可以在使用脚本语言制作算法原型时做到这一点。

因此导数函数最好有签名

void f( double * dy, double * y, double t);

其中维度和常量是外部常量(或将参数数组添加到参数列表)

然后在 rk4 循环中使用它作为

f(k1,  y, t      );
for(i=0; i<n; i++) yt[i] = y[i] + 0.5*h*k1[i];
f(k2, yt, t+0.5*h);
for(i=0; i<n; i++) yt[i] = y[i] + 0.5*h*k2[i];
f(k3, yt, t+0.5*h);
for(i=0; i<n; i++) yt[i] = y[i] +     h*k3[i];
f(k4, yt, t+h);
for(i=0; i<n; i++) y[i] = y[i] + h/6*(k1[i]+2*(k2[i]+k3[i])+k4[i]);

当然,您需要确保每个集成也只分配一次 k1,k2,k3,k4,yt。这可以通过 3 种方式发生:

  • 积分程序包含完整循环,return必要时列出点,
  • 有一些结构包含在第一个集成步骤之前构建的工作数组,
  • 或者您使用较新的动态堆栈分配,它也不应该有过高的分配惩罚。