小时

正式用户

最新动态 3周前

  1. 2月前
    2018-12-07 10:07:54

    @那托儿闹海 小时你的网站流量耗尽不能访问了

    谢谢已经恢复了

  2. 2018-11-29 11:49:13

    @那托儿闹海 请参考
    Barakat, R. (1961). Evaluation of the Incomplete Gamma Function of Imaginary Argument by Chebyshev Polynomials. Mathematics of Computation, 15(73), 7-11. doi:10.2307/2003085

    谢谢哈我看看

  3. 2018-11-24 04:12:52

    @那托儿闹海 可以用Chebyshev插值, 插值还是收敛很快的.

    但这种插值需要节点吗?节点又该怎么求?

  4. 3月前
    2018-11-11 07:10:06
    小时 发表了帖子 如何求任意精度的函数值?

    如果有了任意精度的数据类型,再加上任意精度的加减乘除和开根号运算的方法,请问如何快速求一些特殊函数的任意精度值呢?
    我最关心的是复参数的超几何函数 1F1 (mFn 更好),和 whittaker M 函数,但对普通的 sin,cos, exp 等也感兴趣。
    我试了一下泰勒级数和连续分数,但当 abs(x) 很大的时候收敛太慢,而且误差太大(即中间变量需要的精度需要远大于对函数值的精度要求)。
    我想知道像 Mathematica 这样的软件用什么算法来得到任意精度的? 是用 Arithmetic-Geometric Mean 吗?还是不同的函数不同的变量区间有不同的算法?如果用 AGM 的话,如何获得一个函数的 AGM 表达式?

  5. 4月前
    2018-10-15 10:18:39

    @长岛冰茶 这种量级的个数线性同余大概还能用吧(我们计算物理课检验了一下16807随机数的前$10^9$个

    Numerical Recipes 3ed 一上来就声明
    ”Never use a generator principally based on a linear congruential generator (LCG) or a multiplicative linear congruential generator (MLCG). We say more about this below.”

    基本上这算法早就过时了:
    The idea of LCGs goes back to the dawn of computing, and they were widely used in the 1950s and thereafter. The trouble in paradise first began to be noticed in the mid-1960s

    @ShongLee 如果要比较不同算法的话,这个章节还推荐了其他几种比较安全的算法。 GitHub 上面有现成的代码,应该 c++ 和 fortran 都有。

  6. 2018-10-14 00:04:59

    @ShongLee 我突然想起那个什么复杂度(名字记不住),可以用来衡量序列的随机性,主要思想是“能生成该序列的最短图灵机程序的长度就是这个序列的复杂度”(看到定义或许你知道我说的是什么复杂度)。所以就算伪随机序列周期很大,而代码很短,这不能保证序列的随机性好。这个也可以这样理解,代码越短,包含的信息量越少,产生的序列无论多长,该序列的“信息量”也不可能超越原代码的信息量。

    我劝你还是先试试吧。毕竟这书是数值计算的经典名著, 上面的代码也被检验了几十年了。 不要拿太理论的东西吓唬自己。理论上,足够多个猴子足够长的时间随机敲键盘还能打出莎士比亚的著作呢。

  7. 2018-10-13 13:12:54

    @ShongLee 可以参考 numerical recipes 的 random numbers 章节。里面给出的最好的伪随机数生成器只有十多行代码,周期是 3.13e57,这个算法曾经有一个 1000 刀的奖金奖给从统计规律发现生成的随机数并不完全随机的人,结果无人获奖(现在奖励已经取消了)。这章给出的建议是 Never use the built-in generators in the C and C++ languages(估计 fortran 版也是这么说 fortran 的)。 所以我十分怀疑这个算法不能满足你的需求。

    [attachment:5bc17ebe9723f]
    [attachment:5bc17ec4ba5f6]

  8. 5月前
    2018-08-20 11:44:25
    小时 更新于 请教两个 c++ 问题

    @行人一棹天涯 上面那种写法应该可以吧。C++编译器会选择最匹配的函数。

    非常感谢! 原来是因为 Vector<T>& 要比 T& 更加 specialized。

  9. 6月前
    2018-08-15 10:33:55
    小时 更新于 请教两个 c++ 问题

    另外还有一个关于 template 的问题,比如我另外有一个函数 void func(x, x1, x2), 其中 x 是 Vector 类型, x1 或 x2 可以都为 Vector, 也可以其中一个为普通的 int, double, complex 等标量, 于是我这样声明

    template <class T, class T1, class T2> // 我希望 T, T1, T2 为某种标量类型如 int, double, complex 等
    void func(Vector<T>& x, Vector<T1>& x1, Vector<T2>& x2);
    
    template <class T, class T1, class T2>
    void func(Vector<T>& x, Vector<T1>& x1, T2 x2);
    
    template <class T, class T1, class T2>
    void func(Vector<T>& x, T1 x1, Vector<T2>& x2);

    然而这样写貌似是不行的, 因为例如对于第二个 template, 编译器同样也可以认为 T2 是 Vector<> 类。 而我又不想非常 specific 地声明标量,例如

    template <class T, class T1>
    void func(Vector<T>& x, Vector<T1>& x1, int x2);
    
    template <class T, class T1>
    void func(Vector<T>& x, double x1, Vector<T1>& x2);

    等等。 请问怎么破?

  10. 2018-08-15 03:15:33
    小时 发表了帖子 请教两个 c++ 问题

    我写了一个矢量类 Vector, 一个矩阵类 Matrix, 二者都继承 Base 类。 代码大致如下

    template <class T>
    class Base
    {
        int N; // 矩阵元的总个数
        T* p; // 第一个矩阵元的指针, 动态分配内存储存所有矩阵元
        T& operator()(int i) {return p[i]}
        int size() {return N};
        ...
    }
    
    template <class T>
    class Vector : Base<T>
    {
        void resize(int new_Nsize); // 删除旧内存,分配新大小的内存, N = new_Nsize
        ...
    }
    
    template <class T>
    class Matrix : Base<T>
    {
        int Ncols, Nrows; // 矩阵行数和列数
        void resize(int new_Nrows, int new_Ncols); // 删除旧内存,分配新大小的内存, N = new_Nrows*new_Ncols
        ....
    }

    现在我想写一个支持任何矢量或矩阵的函数, 例如

    template <class T, class T1, class T2>
    void plus(Base<T>& v, Base<T1>& v1, Base<T2>& v2)
    {
        for (int i = 0; i < v1.size(); ++i)
            v(i) = v1(i) + v2(i);
    }

    这个函数运行没有问题. 然而, 我想在函数内部检查 v, v1, v2 是否同为 Vector 或同为 Matrix, 且检查 v1, v2 是否尺寸相同(如果是矩阵,要求行数和列数都相同), 然后再用 resize() 将 v 的尺寸变为和 v1, v2 一样. 这些能做到吗? 要如何实现(尽量保持高 performance)?

查看更多