为防止广告,目前nocow只有登录用户能够创建新页面。如要创建页面请先登录/注册(新用户需要等待1个小时才能正常使用该功能)。

Sgu/111

来自NOCOW
< Sgu
跳转到: 导航, 搜索

高精度开方.可以模拟笔算开方或二分+高精度乘法(需要常数很小).

//P*M
#include<cstdio>
#include<cstring>
#include<cmath>
#include<algorithm>
const long long E8 = 100000000;
int E[8], len, t, l, comp;
char ch1, str1[10000];
using namespace std;
struct bign{
  int l;
  long long a[150];
  bign(){
    l = 0;
    a[0] = 0; 
  }
  bign(int x){
    l = 0;
    a[0] = 0;
    while(x > 0){a[l++] = x % E8; x = x / E8;}
  }
  void prntf(){
    printf("%I64d", a[l-1]);
    for(int i = l-2; i >= 0; i--){
      for(int j = 7; j >=0; j--){
        printf("%I64d", (a[i] % (E[j]*10))/ E[j]);
      }
    }
    printf("\n");
  }
  void div2(){
    for(int i = l-1; i >= 1; i--){
      a[i-1] += (a[i] % 2)*E8;
      a[i] /= 2;
    }
    a[0] = a[0] / 2;
    if(a[l-1] == 0){l--;}
  }
};
int com(bign& a, bign& b){
  if(a.l < b.l){return -1;}
  if(a.l > b.l){return 1;}
  for(int i = a.l-1; i >=0; i--){
    if(a.a[i] < b.a[i]){return -1;}
    if(a.a[i] > b.a[i]){return 1;}
  }
  return 0;
}
bign mid(0), ans(0), r(0), one(1), x1(0);
bign r1(0), r2(0), n(0);
void plus(bign& a, bign& b, bign& c){
  c.a[0] = 0;
  c.l = b.l;
  for(int i = 0; i < c.l; i++){
    c.a[i]+=a.a[i]+b.a[i];
    c.a[i+1] = c.a[i] / E8;
    c.a[i] = c.a[i] % E8;
  }
  for(int i = b.l; i < a.l; i++){
    c.a[i]+=a.a[i];
    c.a[i+1] = c.a[i] / E8;
    c.a[i] = c.a[i] % E8;
  }
  c.l = a.l;
  if(c.a[c.l] != 0){c.l++;}
  return;
}
void minus(bign& a, bign& b, bign& c){
  c.a[0] = 0;
  c.l = b.l; 
  for(int i = 0; i < c.l; i++){
    c.a[i]+=a.a[i]-b.a[i];
    if(c.a[i] < 0){c.a[i]+=E8; c.a[i+1] = -1;}else{c.a[i+1] = 0;}
    c.a[i] = c.a[i] % E8;
  }
  for(int i = c.l; i < a.l; i++){
    c.a[i]+=a.a[i];
    if(c.a[i] < 0){c.a[i]+=E8; c.a[i+1] = -1;}else{c.a[i+1] = 0;}
    c.a[i] = c.a[i] % E8;
  }
  c.l = a.l;
  for(int i = c.l-1; i >= 0; i--){
    if(c.a[i] == 0){c.l--;}else{break;}
  }
  return;
}
void multiply(bign& a, bign& b, bign& d){
  memset(d.a, 0, sizeof(d.a));
  d.a[a.l+b.l-1] = 0;
  for(int i = 0; i < b.l; i++){
    for(int j = 0; j < a.l; j++){
      d.a[i+j] += b.a[i]*a.a[j];
      d.a[i+j+1] += d.a[i+j] / E8;
      d.a[i+j] = d.a[i+j] % E8; 
    }
  }
  if(d.a[a.l+b.l-1] != 0){d.l = a.l+b.l;}else{d.l = a.l+b.l-1;}
  return;
}
int main(){
  freopen("sgu111.in", "r", stdin);
  freopen("sgu111.out", "w", stdout);
  E[0] = 1;
  for(int i = 1; i < 8; i++){E[i] = E[i-1]*10;}
  scanf("%s", &str1);
  len = strlen(str1);
  for(int i = 0; i < len/2; i++){t = str1[i]; str1[i] = str1[len-i-1]; str1[len-i-1] = t;}
  n.l = -1;
  for(int i = 0; i < len; i++){
    if(i % 8 == 0){n.a[++n.l] = str1[i]-48;}else{n.a[n.l] += (str1[i]-48)*E[i%8];}
  }
  r = one;
  l = 0;
  while(r.l < n.l / 2 + 2){plus(r, r, r1); r = r1;}
  r1 = r;
  r2 = r;
  r2.l = max(r2.l-3, 0);//卡下界
  r = one;
  n.l++;
  while(minus(r1, r2, x1), com(x1, r) > 0){
    plus(r1, r2, mid);
    mid.div2();
    multiply(mid, mid, x1);
    comp = com(n, x1);
    if(comp == 0){
      ans = mid;
      break;
    }
    if(comp == -1){
      minus(mid, r, r1);
    }else{
      ans = mid;
      r2 = mid;
    }  
  }
  multiply(r1, r1, x1);
  comp = com(n, x1);
  if(comp == 0){
    ans = r1;
  }
  if(comp == -1){
  }else{
    ans = r1;
  }
  ans.prntf();
  fclose(stdin);
  fclose(stdout);
  return 0;
}
个人工具