Square root digital expansion
Problem 80
It is well known that if the square root of a natural number is not an integer, then it is irrational. The decimal expansion of such square roots is infinite without any repeating pattern at all.
The square root of two is 1.41421356237309504880..., and the digital sum of the first one hundred decimal digits is 475.
For the first one hundred natural numbers, find the total of the digital sums of the first one hundred decimal digits for all the irrational square roots.
简单来说,这题是要我们高精度开方。利用C自带的sqrt库函数肯定是搞笑的,题目要求我们求到根号小数点后100位,double不能达到精度要求。(印象中double的精确值可能只到小数点6、 7位?)
像我们处理大数加法和乘法那样的,同样用数组来表示一个精度足够高的小数,再进行开方的估计,是这题的基本思路。
一、二分逼近
√n的值总是在0~n之间,基于此,而y=√x恰恰是单调函数,于是我们可以用n/2作为初始值,平方一下,比较其与n的大小,如果太大,就让新的估计值t=(n/2+0)/2,如果太小,就让t=(n/2+n)/2,以此类推,直到精度符合要求。
坏处是写小数的高精度加法有点不好写,可以考虑拆分成整数部分和小数部分。或者固定小数点(比如令a[100]为整数部分的个位),在所有的小数末位添0。
二、牛顿法
一个相关的数学定理是:如果g是√n的一个估计值,那么(g+n/g)/2是g的一个更好的估计值。
基于此,可以通过不断迭代来达到精度。
坏处跟上面一样,而且不好估计自己是否达到了精度。
三、手算法
这个方法我也不知道叫什么名字了。。。但思路来源是一个曾经看过的一个叫《如何手算开平方》的文章。
找了篇类似的,链接:https://jingyan.baidu.com/article/4f7d5712f663371a20192797.html
原理什么的百度也有,这也不是数学技术贴,实在是介绍不来。
好处是一次操作我们就可以得到一位的精确值,可以知道什么时候结束,而且写起了也比较简单,只用写个高精度乘法就可以。代码贴下:
#include <stdio.h>
#include <math.h>
#include <string.h>
void tianling(int yu[]);
void found(int yu[],int ans[]);
void xiangcheng(int a[],int b,int ans[]); // a * b = ans
int comp(int a[],int b[]); // a<b 1 a>b 0
void xiangjian(int a[],int b[],int ans[]); // a - b = ans
int main() {
int ans[200]={0},yu[200]={0},j=2,sum=0;
for (int i=2;i<100;i++) {
memset(ans,0,200*sizeof(int));
memset(yu,0,200*sizeof(int));
if (i==j*j) {
j++;
continue;
}
ans[0]=(int)sqrt(i);
yu[0]=i-ans[0]*ans[0];
if (yu[0]>9) {
yu[1]=yu[0]/10;
yu[0]%=10;
}
for (int k=1;k<=100;k++) {
tianling(yu);
found(yu,ans);
}
// if (1) {
// printf("%d\n",i);
// for (int j=100;j>=0;j--) {
// printf("%d",ans[j]);
// }
// printf("\n");
// } //for test
for (int i=1;i<=100;i++) {
sum+=ans[i];
}
if (i==2) printf("%d",sum);
}
printf("%d",sum);
}
int comp(int a[],int b[]) {
int taila=200,tailb=200,max;
while (a[--taila]==0);
while (b[--tailb]==0);
max=taila>tailb?taila:tailb;
for (int i=max;i>=0;i--) {
if (a[i]>b[i]) return 0;
if (a[i]<b[i]) return 1;
}
}
void xiangjian(int a[],int b[],int ans[]) {
int taila=200,tailb=200,max;
memset(ans,0,120*sizeof(int));
while (a[--taila]==0);
while (b[--tailb]==0);
max=taila>tailb?taila:tailb;
for (int i=0;i<=taila;i++) {
ans[i]+=a[i]-b[i];
if (ans[i]<0) {
ans[i]+=10;
ans[i+1]--;
}
}
}
void tianling(int yu[]) {
int i=199;
while (i>=0) {
if (yu[i]!=0) break;
i--;
}
for (int j=i;j>=0;j--) {
yu[j+2]=yu[j];
}
yu[0]=yu[1]=0;
}
void xiangcheng(int a[],int b,int ans[]) {
memset(ans,0,200*sizeof(int));
for (int i=0;i<200;i++) {
ans[i]=a[i]*b;
}
for (int i=0;i<200;i++) {
if (ans[i]>9) {
if (i==199) {
printf("overflow");
break;
}
ans[i+1]+=ans[i]/10;
ans[i]%=10;
}
}
}
void found(int yu[],int ans[]) {
int temp[200]={0},newyu[200]={0},newans[200]={0},k;
xiangcheng(ans,20,newans);
for (k=9;k>=0;k--) {
newans[0] = k;
memset(temp,0,200*sizeof(int));
xiangcheng(newans,k,temp);
if (comp(yu,temp)==0) {
break;
}
}
xiangjian(yu,temp,newyu);
for (int i=0;i<200;i++) {
yu[i] = newyu[i];
}
xiangcheng(ans,10,newans);
newans[0] = k;
for (int i=0;i<200;i++) {
ans[i] = newans[i];
}
}