多参数和谐波函数如何查询GSL ODE
How to consult GSL ODE with many parameters and harmonic functions
我正在使用 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;
}
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 使用前。
如果您使用的是参数,它们很可能属于同一类型 (double)。在那种情况下,这也可以使用数组来解决,然后从 func and/or jac.
访问元素
另一种选择是使用 gsl_vector,然后“获取”函数内的值。这将涉及免费使用。
我正在使用 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;
}
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 使用前。
如果您使用的是参数,它们很可能属于同一类型 (double)。在那种情况下,这也可以使用数组来解决,然后从 func and/or jac.
访问元素另一种选择是使用 gsl_vector,然后“获取”函数内的值。这将涉及免费使用。