开心一刻，考考你的计算能力

G-Spider

7# mathe 我用matlab算了一下，最后一行数据不一样。 --------------------------------------------------------------------- 。。。 3799248216530963946104908731687 718950889

sheng_jianguo

7# mathe 我用matlab算了一下，最后一行数据不一样。 --------------------------------------------------------------------- 。。。 3799248216530963946104908731687 718950889 G-Spider 发表于 2011-1-13 19:30

用matlab算，怎样才能得到如此高的精度？（我也用matlab算，好像精度不够）

G-Spider

12# sheng_jianguo >>help vpa VPA Variable precision arithmetic. R = VPA(S) numerically evaluates each element of the double matrix S using variable precision floating point arithmetic with D decimal digit accuracy, where D is the current setting of DIGITS. The resulting R is a SYM. VPA(S,D) uses D digits, instead of the current setting of DIGITS. D is an integer or the SYM representation of a number. It is important to avoid the evaluation of an expression using double precision floating point arithmetic before it is passed to VPA. For example, phi = vpa((1+sqrt(5))/2) first computes a 16-digit approximation to the golden ratio, then converts that approximation to one with d digits, where d is the current setting of DIGITS. To get full precision, use unevaluated string or symbolic arguments, phi = vpa('(1+sqrt(5))/2') or s = sym('sqrt(5)') phi = vpa((1+s)/2); ////////////////////////// 可见matlab有良好的大数运算功能。 vpa('173746*sin(10^22)+94228*log(17.1)-78487*exp(0.42)',1000) %注意要加上引号，不然每一项不到15位的精度，最后的结果比1000的精度略低一点，存在舍入

sheng_jianguo

13# G-Spider 谢谢！

winxos

我是打酱油的， 过来看看各位老大。

wayne

9# wayne 代码不好，在windows的MinGW下编译，即使开了-Wall也不警告，运行也很正常。 在linux下编译没问题，但运行时会出现segmentation fault . 检查了很长时间，原来是argv[1] 使然。 现在更正如下：

1. #include <stdio.h>
2. #include <stdlib.h>
3. #include <mpfr.h>
4. int main (int argc, char *argv[])
5. {int m;
6. if (argc>1) m= atoi(argv[1])-1;
7. else printf("usage: %s n\n",argv[0]),exit(1);
8. mp_prec_t p = 100+5*m;
9. mpfr_t a, b, c, d;
10. mpfr_inits2 (p, a, b, c, d, (mpfr_ptr) 0);
11. mpfr_set_str (a, "1e22", 10, GMP_RNDN);
12. mpfr_sin (a, a, GMP_RNDN);
13. mpfr_mul_ui (a, a, 173746, GMP_RNDN);
14. mpfr_set_str (b, "17.1", 10, GMP_RNDN);
15. mpfr_log (b, b, GMP_RNDN);
16. mpfr_mul_ui (b, b, 94228, GMP_RNDN);
17. mpfr_set_str (c, "0.42", 10, GMP_RNDN);
18. mpfr_exp (c, c, GMP_RNDN);
19. mpfr_mul_si (c, c, -78487, GMP_RNDN);
20. mpfr_add (d, a, b, GMP_RNDN);
21. mpfr_add (d, d, c, GMP_RNDN);
22. mpfr_printf ("d=%.*Re\n",m,d);
23. mpfr_clears (a, b, c, d, NULL);
24. return 0;
25. }

复制代码