标签: numerical-analysis

我正在使用BigDecimal进行简单的乘法运算,并且在乘以零时我发现了一些奇怪的行为(在这个用例中乘以零是正确的).

基本数学告诉我,任何乘以零的数字都将等于零(参见:零产品属性和乘法属性)

但是,以下代码将始终失败并出现相同的错误:

assertEquals(new BigDecimal(0), new BigDecimal(22.3).multiply(new BigDecimal(0)));

java.lang.AssertionError: 
Expected :0
Actual   :0E-48

Run Code Online (Sandbox Code Playgroud)

这是BigDecimal的不准确,还是我在某处遗漏了数学的一些利基分支？

注意:在IntelliJ 11中运行的JDK 1.6.0_27

我想计算时间序列的功率谱密度; 做一些带通,低通和高通滤波; 也许其他一些基本的东西.

有没有一个很好的开源Java库来做到这一点？

我有点没有成功(例如,谷歌搜索"功率谱密度java"或"信号处理java"和点击链接,查看Apache Commons,Sourceforge,java.net等).

有许多小程序,书籍,教程,商业产品等,不符合我的需求.

更新:我找到了 org.apache.commons.math.transform进行傅立叶变换.这不实现功率谱密度,带通等,但它是一些东西.

我很想知道R的平均函数使用什么算法.是否有一些参考该算法的数值属性？

我在summary.c中找到了以下C代码:do_summary():

case REALSXP:
PROTECT(ans = allocVector(REALSXP, 1));
for (i = 0; i < n; i++) s += REAL(x)[i];
s /= n;
if(R_FINITE((double)s)) {
    for (i = 0; i < n; i++) t += (REAL(x)[i] - s);
    s += t/n;
}
REAL(ans)[0] = s;
break;

它似乎做了一个直线意思:

for (i = 0; i < n; i++) s += REAL(x)[i];
s /= n;

然后它添加了我假设的数值修正,它似乎是与数据平均值的平均差异:

for (i = 0; i < n; i++) t += (REAL(x)[i] - s);
s += t/n;

我无法在任何地方跟踪此算法(意味着不是一个很棒的搜索词).

任何帮助将非常感激.

我正在寻找64位(无符号)立方根的快速代码.(我正在使用C并使用gcc进行编译,但我认为所需的大部分工作都是语言和编译器无关的.)我将通过ulong表示64位未分配的整数.

给定输入n我需要(整数)返回值r

r * r * r <= n && n < (r + 1) * (r + 1) * (r + 1)

也就是说,我想要n的立方根,向下舍入.基本代码如

return (ulong)pow(n, 1.0/3);

因为向范围的末端舍入而不正确.简单的代码就像

ulong
cuberoot(ulong n)
{
    ulong ret = pow(n + 0.5, 1.0/3);
    if (n < 100000000000001ULL)
        return ret;
    if (n >= 18446724184312856125ULL)
        return 2642245ULL;
    if (ret * ret * ret > n) {
        ret--;
        while (ret * ret * ret > n)
            ret--;
        return ret;
    }
    while ((ret + 1) * (ret + …

我在几个地方(在文档和这篇博文中阅读:http://blogs.mathworks.com/loren/2007/05/16/purpose-of-inv/)中读到,在Matlab中使用inv不是建议因为它很慢而且不准确.

我试图找出这种不准确的原因.截至目前,谷歌没有给出有趣的结果,所以我认为有人可以指导我.

谢谢 !

对于数学分析的大学课程,我们正在从Maple过渡到Numpy和Sympy的组合,以获得课程材料的各种插图.这是因为学生们之前已经学习过Python.

我们遇到的困难之一是在Python中模拟固定精度.Maple允许用户指定小数精度(比如10或20位),然后在每次计算时都使用该精度,这样您就可以看到舍入误差的影响.在Python中,我们尝试了一些方法来实现这一点:

• Sympy具有指定位数的舍入函数.
• Mpmath支持自定义精度.

然而,这不是我们想要的.这些选项计算确切结果并将精确结果四舍五入到指定的位数.我们正在寻找一种能够以指定精度进行每次中间计算的解决方案.例如,可以显示分割两个非常小的数字时可能发生的舍入误差的东西.

到目前为止,最好的解决方案似乎是Numpy中的自定义数据类型.使用float16,float32和float64,我们至少能够指出可能出现的问题.这里的问题是我们总是需要使用一个元素的数组,并且我们仅限于这三种数据类型.

为了我们的目的,还有更灵活的东西吗？或者我们正在寻找隐藏在mpmath文档中某处的东西？当然,通过在舍入函数中包装计算的每个元素有一些变通方法,但这会使代码模糊到学生.

我正在尝试使用牛顿方法的不同实现来计算平方根.一个重要的决定是什么时候终止算法.

显然,它不会做使用之差的绝对值y*y和x,其中y是的平方根的当前估计x,因为对于大的值x也可能无法代表具有足够精度的平方根.

所以我应该使用相对标准.天真的我会用这样的东西:

static int sqrt_good_enough(float x, float y) {
  return fabsf(y*y - x) / x < EPS;
}

这看起来效果很好.但是最近我开始阅读Kernighan和Plauger的编程风格元素,他们在第1章中给出了相同算法的Fortran程序,其终止标准(用C翻译)将是:

static int sqrt_good_enough(float x, float y) {
  return fabsf(x/y - y) < EPS * y;
}

两者在数学上是等价的,但是有理由选择一种形式而不是另一种吗？

我注意到如果A是NxN矩阵并且它具有逆矩阵.但是inv()和pinv()函数输出是不同的. - 我的环境是Win7x64 SP1,Matlab R2012a,Cygwin Octave 3.6.4,FreeMat 4.2

看看Octave的例子:

A = rand(3,3)
A =
0.185987   0.192125   0.046346
0.140710   0.351007   0.236889
0.155899   0.107302   0.300623

pinv(A) == inv(A)
ans =
0 0 0
0 0 0
0 0 0

• 这一切都同ans通过上面运行在Matlab同一命令的结果.

• 我计算inv(A)*A或A*inv(A),结果是Octave和Matlab中的身份3x3矩阵.
• 结果A*pinv(A)和pinv(A)*A在Matlab和FreeMat身份3x3矩阵.
• 的结果 A*pinv(A)是Octave中的身份3x3矩阵.
• 结果pinv(A)*A是不是在倍频身份3x3矩阵.

我不知道原因 inv(A) != pinv(A),我已经考虑了矩阵中元素的细节.这似乎是导致这个问题的浮动精度问题.

点点后10位数可能不同,如下所示:

• 6.65858991579923298331777914427220821380615200000000inv(A)(1,1)反对的元素

• 6.65858991579923209513935944414697587490081800000000 元素 pinv(A)(1,1)

http://en.wikipedia.org/wiki/Binary_GCD_algorithm

此维基百科条目有一个非常不满意的含义:二进制GCD算法是在同一时间高达60%,比标准的Euclid算法更有效,但迟至1998年克努特的结论是,只有在效率15%的涨幅对他的当代电脑.

还有15年过去了......今天这两种算法如何与硬件的进步相结合？

二进制GCD在低级语言中是否继续优于欧几里得算法,但由于其在Java等高级语言中的复杂性而萎缩？或者是现代计算的差异？

我为什么关心你可能会问？我今天碰巧需要像今天的1000亿这样处理:)这是一个生活在计算时代(糟糕的欧几里德)的干杯.

我正在研究一个模拟系统.我将很快获得实验数据(直方图),用于几个模拟输入的实际值分布.

当模拟运行时,我希望能够生成与测量分布匹配的随机值.我宁愿这样做而不存储原始直方图.有什么好方法

1. 将直方图映射到表示分布的一组参数？
2. 在运行时生成基于这些参数的值？

编辑:输入数据是几种不同类型事件的事件持续时间.我希望不同的类型具有不同的分布函数.