#include <iostream>
#include <cmath>
using namespace std;
int main()
{
double a = sqrt(2);
cout << a << endl;
}
的許多數字你好,這是找開方的2它打印只是1.41421輸出如何實現它的方式,使得其將打印小數點後20萬位程序點發現2的平方根儘可能
1.41421 ..........高達200 000位
是否有任何的方法來打印這樣嗎?
#include <iostream>
#include <cmath>
using namespace std;
int main()
{
double a = sqrt(2);
cout << a << endl;
}
的許多數字你好,這是找開方的2它打印只是1.41421輸出如何實現它的方式,使得其將打印小數點後20萬位程序點發現2的平方根儘可能
1.41421 ..........高達200 000位
是否有任何的方法來打印這樣嗎?
這裏是你的問題使用GNU GMP庫的代碼。該代碼將打印1000位數的sqrt(2),增加註釋行中的數字以滿足您的請求。
#include <stdio.h>
#include <gmp.h>
int main(int argc, char *argv[])
{
mpf_t res,two;
mpf_set_default_prec(1000000); // Increase this number.
mpf_init(res);
mpf_init(two);
mpf_set_str(two, "2", 10);
mpf_sqrt (res, two);
gmp_printf("%.1000Ff\n\n", res); // increase this number.
return 0;
}
請用下面的命令進行編譯:
$gcc gmp.c -lgmp -lm -O0 -g3
就雙精度算法的精度而言,您給出的例子是準確的,這是大多數C++編譯器使用的最高精度。一般而言,計算機不適合做更高精度的計算。 如果這是某種家庭作業,那麼我懷疑你需要找出一個計算算法 - 你需要保持你自己的數字數組,以保持你所需要的所有精度。 如果您有一些真實世界的應用程序,您絕對應該使用專門用於執行此類算術的高精度庫(GMP是一個很好的開源可能性) - 這是一個複雜的輪子,不需要重新創建。
sqrt(2) = (239/169)*1/sqrt(1-1/57122)
和1/SQRT(1-1/57122)可以有效地利用泰勒級數來計算擴展:
1/sqrt(1-x) = 1 + (1/2)x + (1.3)/(2.4)x^2 + (1.3.5)/(2.4.6)x^3 + ...
There's also a C program available that uses this method(我稍微重新格式化並糾正它):
/*
** Pascal Sebah : July 1999
**
** Subject:
**
** A very easy program to compute sqrt(2) with many digits.
** No optimisations, no tricks, just a basic program to learn how
** to compute in multiprecision.
**
** Formula:
**
** sqrt(2) = (239/169)*1/sqrt(1-1/57122)
**
** Data:
**
** A big real (or multiprecision real) is defined in base B as:
** X = x(0) + x(1)/B^1 + ... + x(n-1)/B^(n-1)
** where 0<=x(i)<B
**
** Results: (PentiumII, 450Mhz)
**
** 1000 decimals : 0.02seconds
** 10000 decimals : 1.7s
** 100000 decimals : 176.0s
**
** With a little work it's possible to reduce those computation
** times by a factor of 3 and more.
*/
#include <stdio.h>
#include <stdlib.h>
long B = 10000; /* Working base */
long LB = 4; /* Log10(base) */
/*
** Set the big real x to the small integer Integer
*/
void SetToInteger(long n, long* x, long Integer)
{
long i;
for (i = 1; i < n; i++)
x[i] = 0;
x[0] = Integer;
}
/*
** Is the big real x equal to zero ?
*/
long IsZero(long n, long* x)
{
long i;
for (i = 0; i < n; i++)
if (x[i])
return 0;
return 1;
}
/*
** Addition of big reals : x += y
** Like school addition with carry management
*/
void Add(long n, long* x, long* y)
{
long carry = 0, i;
for (i = n - 1; i >= 0; i--)
{
x[i] += y[i] + carry;
if (x[i] < B)
carry = 0;
else
{
carry = 1;
x[i] -= B;
}
}
}
/*
** Multiplication of the big real x by the integer q
*/
void Mul(long n, long* x, long q)
{
long carry = 0, xi, i;
for (i = n - 1; i >= 0; i--)
{
xi = x[i] * q;
xi += carry;
if (xi >= B)
{
carry = xi/B;
xi -= carry * B;
}
else
carry = 0;
x[i] = xi;
}
}
/*
** Division of the big real x by the integer d
** Like school division with carry management
*/
void Div(long n, long* x, long d)
{
long carry = 0, xi, q, i;
for (i = 0; i < n; i++)
{
xi = x[i] + carry * B;
q = xi/d;
carry = xi - q * d;
x[i] = q;
}
}
/*
** Print the big real x
*/
void Print(long n, long* x)
{
long i;
printf("%ld.", x[0]);
for (i = 1; i < n; i++)
printf("%04ld", x[i]);
printf("\n");
}
/*
** Computation of the constant sqrt(2)
*/
int main(void)
{
long NbDigits = 200000, size = 1 + NbDigits/LB;
long* r2 = malloc(size * sizeof(long));
long* uk = malloc(size * sizeof(long));
long k = 1;
/*
** Formula used:
** sqrt(2) = (239/169)*1/sqrt(1-1/57122)
** and
** 1/sqrt(1-x) = 1+(1/2)x+(1.3)/(2.4)x^2+(1.3.5)/(2.4.6)x^3+...
*/
SetToInteger(size, r2, 1); /* r2 = 1 */
SetToInteger(size, uk, 1); /* uk = 1 */
while (!IsZero(size, uk))
{
Div(size, uk, 57122); /* uk = u(k-1)/57122 * (2k-1)/(2k) */
Div(size, uk, 2 * k);
Mul(size, uk, 2 * k - 1);
Add(size, r2, uk); /* r2 = r2+uk */
k++;
}
Mul(size, r2, 239);
Div(size, r2, 169); /* r2 = (239/169)*r2 */
Print(size, r2); /* Print out of sqrt(2) */
free(r2);
free(uk);
return 0;
}
大約需要一分鐘來計算開方200,000位(2)。
但是,請注意,由於累計舍入錯誤,所產生的最後11位數字不正確,如果需要200,000位正確數字,則需要運行200,012位數字。
這是一個解決方案,它能夠在不到一分鐘的時間內以良好的舊版Prolog編程語言計算sqrt(2)的一百萬位數。它是基於求解Pell方程,也參見here:
p*p+1 = 2*q*q
的重複週期關係P'= 3P + 4Q和q'= 2P + 3Q可以鑄造成矩陣乘法。即我們發現,如果我們乘以向量[P,Q]由係數矩陣我們得到的矢量[P「 Q」]:
| p' | | 3 4 | | p |
| | = | | * | |
| q' | | 2 3 | | q |
對於矩陣A,我們可以使用一個二進制遞歸使得我們可以計算O(log n)操作中的A^n。我們將需要大量的數量。我用這實驗代碼here由此主程序然後簡單地:
/**
* pell(N, R):
* Compute the N-the solution to p*p+1=2*q*q via matrices and return
* p/q in decimal format in R.
*/
pell(N, R) :-
X is [[3,4],[2,3]],
Y is X^N,
Z is Y*[[1],[1]],
R is Z[1,1]*10^N//Z[2,1].
下面的屏幕截圖示出了定時和一些結果。我使用了10次100萬次迭代。可以將結果與本頁面here進行比較。
缺少的請告訴我仍然是一個明確的標準和計算,說多的數字是如何保持穩定。我們需要更多時間來做到這一點。
編輯2016年12月20日:
我們通過一個上限的相對誤差的改進的代碼一點點,並通過輕推結果進一步計算穩定數字。 100萬位的計算時間現在低於2秒:
?- time(pell(653124, _, S)).
% Uptime 1,646 ms, GC Time 30 ms, Thread Cpu Time 1,640 ms
S = -1000000
'double'會給你大約16位有效數字。如果你想要200000,你需要一個任意的精確庫(例如GMP)。 – 2013-03-12 13:08:07
http://en.wikipedia.org/wiki/Methods_of_computing_square_roots – 2013-03-12 13:08:30
我不認爲使用內置的sqrt函數會做你想做的。我猜這是一項家庭作業。 – crush 2013-03-12 13:11:54