多参数和谐波函数如何查询GSL ODE

Question

我正在使用 GSL 研究非线性微分方程。问题是我对 C 的东西很陌生。我刚刚将 GNU 站点上的样本改编成了我现在感兴趣的方程式。

这是等式：

d2x/dt2 + r*dx/dy + cos(x) + v*cos(2*x+0.4) E1*sin(wt) + E2*sin(2*w*t+a) = 0

我卡住的是我不知道如何在代码中插入多个参数。此外，我不知道如何在这段代码中使用余弦或正弦函数。

我一直在 Google 上搜索，试图找出这个问题。我找不到任何对我有帮助的东西。

#include <stdio.h>
#include <gsl/gsl_errno.h>
#include <math.h>
#include <gsl/gsl_math.h>
#include <gsl/gsl_matrix.h>
#include <gsl/gsl_odeiv2.h>

int func (double t, const double x[], double y[], void *params)
{
    double r = *(double *)params;
    double v = *(double *)params;
    double w = *(double *)params;
    double E1 = *(double *)params;
    double E2  = *(double *)params;
    double a  = *(double *)params;
    y[0] = x[1];
    y[1] = -r*x[1] - cos(x[0]) - v*cos(2*x[0]+0.4) - E1*sin(w*t) - E2*sin(2*w*t+a);
    return GSL_SUCCESS;
}

int jac (double t, const double x[], double *dydx, double dydt[], void *params)
{
    double r = *(double *)params;
    double v = *(double *)params;
    double w = *(double *)params;
    double E1 = *(double *)params;
    double E2  = *(double *)params;
    double a  = *(double *)params;
    gsl_matrix_view dydx_mat = gsl_matrix_view_array (dydx, 2, 2);
    gsl_matrix * m = &dydx_mat.matrix;
    gsl_matrix_set (m, 0, 0, 0.0);
    gsl_matrix_set (m, 0, 1, 1.0);
    gsl_matrix_set (m, 1, 0, sin(x[0]) + 2*v*sin(2*x[0]+0.4));
    gsl_matrix_set (m, 1, 1, -r);
    dydt[0] = 0.0;
    dydt[1] = 0.0;
    return GSL_SUCCESS;
}

int main (void)
{
    double r = 0.0;
    double v = 0.0;
    double w = 2.4;
    double E1 = -2.3;
    double E2 = 0;
    double a = 0.7;
    gsl_odeiv2_system sys = {func, jac, 2, &r, &v, &w, &E1, &E2, &a};

    gsl_odeiv2_driver *d = gsl_odeiv2_driver_alloc_x_new (&sys, gsl_odeiv2_step_rk8pd, 1e-6, 1e-6, 0.0);
    int i;
    double t = 0.0, t1 = 10000;
    double x[2] = {0.0, 0.0};

    for (i = 1 ; i<=10000; i++)
        {
            double ti = i*t1/10000;
            int status = gsl_odeiv2_driver_apply (d, &t, ti, x);

            if (status != GSL_SUCCESS)
                {
                    printf("error, return value%d\n", status);
                    break;
                }
            printf("%.5e %.5e %.5e\n", t, x[0], x[1]);
        }

    gsl_odeiv2_driver_free (d);
    return 0;
}

Answer 1

params 参数是指向某些任意数据结构 的指针（地址/内存位置）。在 GSL 文档的示例中，他们的等式仅包含一个参数，这意味着可以只传递 double 精度数字的地址。

但是，对于您的问题，您需要访问 6 个不同的参数。您不能访问具有相同地址的每个参数！

/* this doesn't work! */
double r  = *(double *)params;
double v  = *(double *)params;
double w  = *(double *)params;
double E1 = *(double *)params;
double E2 = *(double *)params;
double a  = *(double *)params;

由于所有地址都相同，所以您指的是同一个号码。要解决此问题，您可以：将所有参数存储在长度为 6 的数组中，或者将它们存储在预定义的数据结构中。后一种方法更具可读性，所以我将证明这一点。

首先定义一个数据类型来指定要存储的参数：

struct param_type {
    double r;
    double v;
    double w;
    double E1;
    double E2;
    double a;
};

现在，在main函数中创建一个这种类型的结构并存储参数的实际值：

struct param_type my_params = {r, v, w, E1, E2, a};

定义系统时，您存储一个指向该系统的指针 struct param_type:

gsl_odeiv2_system sys = {func, jac, 2, &my_params};

要在 func 和 jac 中使用参数，您只需将 params 参数从通用指针 (void *) 转换为特定数据的指针输入 (struct param_type *):

struct param_type *my_params_pointer = params;

（请注意，在 C++ 中，这必须使用 显式转换 。）最后，您可以通过以下方式访问参数：

double r  = my_params_pointer->r;
double v  = my_params_pointer->v;
double w  = my_params_pointer->w;
double E1 = my_params_pointer->E1;
double E2 = my_params_pointer->E2;
double a  = my_params_pointer->a;

这里用箭头->代替点.因为my_params_pointer是一个指针需要dereferenced 使用前。

Answer 2

如果您使用的是参数，它们很可能属于同一类型 (double)。在那种情况下，这也可以使用数组来解决，然后从 func and/or jac.

访问元素

另一种选择是使用 gsl_vector，然后“获取”函数内的值。这将涉及免费使用。

多参数和谐波函数如何查询GSL ODE

How to consult GSL ODE with many parameters and harmonic functions

c

gsl

differential-equations

nonlinear-functions