Tag: gsl

使用gsl编写Runge-Kutta ODE求解器

自从我做了任何C / c ++以来已经有一段时间了，但我想用gsl库编写一个ODE求解器来解决下面的ODE集 $$ u'(r)=up(r)$$ $$ up'(r)=-(2*(r-1)/(r*(r-2)))*up(r)-((r*r/((r-2)*(r-2)))-(2/r*(r-2)))*u(r) $$ 所以在gsl表示法中我的y [0] = u，y [1] == up，并且上面的RHS定义了f [0]和f [1]。 然后可以从这些定义中计算雅可比行列式和dfdr（通常它们的“时间”变量称为“t”而不是“r”）。 这样做的原因是因为我对Mathematica有速度问题。 我在ODE解算器的文档末尾使用了gsl示例代码，并尝试将其调整到我的问题，如下所示： #include #include #include #include #include int func (double r, const double y[], double f[], void *params) { double mu = *(double *)params; f[0] = y[1]; f[1] = -(2*(r-1)/(r*(r-2)))*y[1]-((r*r/((r-2)*(r-2)))-(2/r*(r-2)))*y[0]; return GSL_SUCCESS; } /* void tester […]

使用64位GCC在Cygwin上编译64位GSL

我试图在Windows 7机器上编译64位GSL。 以下是我采取的步骤： 下载并解压缩此处的GSL 1.15源代码。 通过编译最小程序测试我在Cygwin shell中有64位版本的GCC // simple.C int main() { ; return 0; } 运用 x86_64-w64-mingw32-gcc -m64 simple.C -o simple 在untarred文件夹中，我想将x86_64-w64-mingw32-gcc编译器传递给./configure但我不确定如何。 我查看了configure文件，但这很大，似乎是使用autoconf生成的。

使用GSL计算矩阵Kronecker乘积的有效方法

我的算法的瓶颈是我的functionKronecker产品叫KPro： gsl_matrix *KPro(gsl_matrix *a, gsl_matrix *b) { int i, j, k, l; int m, p, n, q; m = a->size1; p = a->size2; n = b->size1; q = b->size2; gsl_matrix *c = gsl_matrix_alloc(m*n, p*q); double da, db; for (i = 0; i < m; i++) { for (j = 0; j < p; j++) { […]

gsl复杂矩阵*复杂向量

有人可以帮我弄清楚如何进行这项操作。 我正在查看文档，但它已经很久了，因为我必须做任何线性代数类型的东西，我有点迷失。 我有一个叫做Y的14×14复杂矩阵和一个叫做I的复矢量。我需要将它们相乘并将结果设置为一个名为IL的复矢量。 到目前为止，我已经发现我需要使用： gsl_blas_zgemv (CBLAS_TRANSPOSE_t TransA, const gsl_complex alpha, const gsl_matrix_complex * A, const gsl_vector_complex * x, const gsl_complex beta, gsl_vector_complex * y) 但我不确定究竟发生了什么。 不知道到底发生了什么。 像这样的东西？ 但是alpha和beta是什么？ gsl_blas_zgemv(CblasNoTrans, ???, &Y, &I, ???, IL);

OpenMP和GSL RNG – 性能问题 – 4个线程实现比纯序列1（四核CPU）慢10倍

我正在尝试将我的C项目从顺序编程转换为并行编程。 尽管大多数代码现在已经从头开始重新设计，但随机数的生成仍然是其核心。 因此，随机数发生器（RNG）的不良性能会严重影响程序的整体性能。 我写了一些代码行（见下文），以显示我面临的问题而没有太多冗长。 问题如下：每次线程数增加时，性能都会明显变差。 在这个工作站（linux内核2.6.33.4; gcc 4.4.4; intel四核CPU）中，无论迭代次数n多少，并行for循环使用nt = 4比使用nt = 1大约长10倍。 这种情况似乎在这里描述，但焦点主要集中在fortran，一种我不太了解的语言，所以我非常感谢一些帮助。 我试图按照他们的想法创建不同的RNG（使用不同的种子）来访问每个线程，但性能仍然很差。 实际上，每个线程的这个不同的播种点也让我感到困惑，因为我无法看到最终如何保证生成的数字的质量（缺乏相关性等）。 我已经考虑过完全放弃GSL并自己实现一个随机生成器算法（例如Mersenne-Twister），但我怀疑我稍后会遇到同样的问题。 非常感谢您的回答和建议。 请问我可能忘记提及的任何重要事项。 编辑：由lucas1024（pragma for-loop声明）和JonathanDursi（播种;将“a”设置为私有变量）建议的更正。 multithreading模式下的性能仍然非常低迷。 编辑2： Jonathan Dursi建议的实施解决方案（见评论）。 #include #include #include #include #include double d_t (struct timespec t1, struct timespec t2){ return (t2.tv_sec-t1.tv_sec)+(double)(t2.tv_nsec-t1.tv_nsec)/1000000000.0; } int main (int argc, char *argv[]){ double a, b; int i,j,k; int n=atoi(argv[1]), seed=atoi(argv[2]), […]

浮点优化 – 指南

我们需要通过在C / C ++中实现特定算法来解决的大多数科学计算问题要求精度远低于双精度。 例如， 1e-6 1e-7精度覆盖了99%的ODE求解器或数值积分的情况。 即使在我们确实需要更高精度的极少数情况下，通常数值方法本身也会失败，然后才能达到接近双精度的精度。 示例：即使在因舍入误差而求解标准nostiff常微分方程时，我们也不能期望从简单的Runge-Kutta方法获得1e-16精度。 在这种情况下，双精度要求类似于要求更好地逼近错误答案。 然后，在大多数情况下，积极的浮点优化似乎是一个双赢的局面，因为它使您的代码更快（更快！）并且它不会影响您的特定问题的目标准确性。 也就是说，确保特定的实现/代码对fp优化是稳定的似乎非常困难。 古典（有点令人不安）的例子：GNL，GNU科学图书馆，不仅是市场上的标准数字图书馆，而且它是一个写得很好的图书馆（我无法想象自己做得更好）。 但是，GSL对fp优化不稳定。 实际上，例如，如果使用intel编译器编译GSL，那么除非打开-fp-model strict标志关闭fp优化，否则它的内部测试将失败。 因此，我的问题是：是否存在编写代码的一般准则，这些准则对于积极的浮点优化是稳定的。 这些指南语言（编译器）是否具体。 如果是这样，那么C / C ++（gcc / icc）的最佳实践是什么？ 注1：这个问题不是询问gcc / icc中的fp优化标志是什么。 注2：这个问题并不是要求关于C / C ++优化的一般指导原则（比如不要对很多小型函数使用虚函数）。 注3：这个问题并不是要求大多数标准fp优化列表（比如x / x – > 1）。 注4：我坚信这不是类似于传统的“最酷的服务器名称”的主观/偏离主题的问题。 如果您不同意（因为我没有提供具体示例/代码/问题），请将其标记为社区维基。 我对答案比对获得一些状态点更感兴趣（不是它们不重要 – 你明白了！）。

使用GPU随机数

我正在研究使用nvidia GPU进行蒙特卡罗模拟。 但是，我想使用gsl随机数生成器以及并行随机数生成器，如SPRNG。 有谁知道这是否可能？ 更新 我使用GPU玩过RNG。 目前还没有一个很好的解决方案。 SDK附带的Mersenne Twister并不适合（我的）Monte-Carlo模拟，因为生成种子需要相当长的时间。 NAG库更有前途。 您可以批量生成RN，也可以在单个线程中生成RN。 但是，目前仅支持少数分布 – Uniform，exponential和Normal。

使用GSL的双摆解决方案

我正在学习使用GSL来解决ODE问题。 我想用GSL ODE函数解决双摆问题。 我决定使用这个等式： （来源： http ： //www.physics.usyd.edu.au/~wheat/dpend_html/ ） 我的代码： #include #include #include “gsl/gsl_errno.h” #include “gsl/gsl_matrix.h” #include “gsl/gsl_odeiv2.h” #include “constants.h” double L1; double L2; double M1; double M2; double T_START; double T_END; double S1_ANGLE; double S2_ANGLE; double V1_INIT; double V2_INIT; int func(double t, const double y[], double f[], void *params) { /* * y[0] = […]