[手動][大數][三角函數] sin x 運算＠Morris' Blog

2011-08-02 21:50:36| 人氣895| 回應0 | 上一篇 | 下一篇

[手動][大數][三角函數] sin x 運算

推薦 0 收藏 0 轉貼0 訂閱站台

sin x = x - x³ / 3! + x⁵ / 5! - x⁷ / 7! + ... for all x.
x 單位是弧度,
下面程式碼是我自己做的, 有錯誤的話, 請多多包涵, 由於只會做到 x^499/499!
所以 x 不能太大, 不然收斂不到, 我沒有做同餘的處理, 請多多包涵
範例輸入:
-0.523593766666666666666666666666666666666666666666666666 //-30度
50
範例輸出:
-0.49999566213168835891110783336543526239570962009857942732098781699100110611002
57315819223394605303886564414913897841757255904112772246821605497279605025577762
9952351298307494694247274622940860856025771
-0.26237485370392878591439364691262254588666596497123852294393559241445895258995
07517889236912853523696585990137026532086111068236521550086251743962915519539144
6454176291446477811029972808232469547290065

#include<stdio.h>
#include<stdlib.h>
#include<string.h>
#include<math.h>
#define L 1000
typedef struct {
   char neg;
   int pv[L], pl;
   int hv[L], hl;
}BigNumber;
/*
123.456
pv[] = '321'
hv[] = '456'
*/
BigNumber minus(BigNumber , BigNumber );
int compare(BigNumber A, BigNumber B) {/*兩個都是正數*/
   int a, h;
   if(A.pl > B.pl) return 1;
   if(A.pl < B.pl) return 0;
   for(a = A.pl; a >=0; a--)
       if(A.pv[a] > B.pv[a])
           return 1;
       else if(A.pv[a] < B.pv[a])
           return 0;
   for(a = 0; a <= A.hl && a <= B.hl; a++)
       if(A.hv[a] > B.hv[a])
           return 1;
       else if(A.hv[a] < B.hv[a])
           return 0;
   if(A.hl > B.hl) return 1;
   else if(A.hl < B.hl) return 0;
   return 1;
}
BigNumber addit(BigNumber A, BigNumber B) {
   if(A.neg == 0 && B.neg == 1) {
       B.neg = 0;return minus(A, B);
   }
   else if(A.neg == 1 && B.neg == 0) {
       A.neg = 0;return minus(B, A);
   }
   int a, p, h;
   BigNumber C;
   C.neg = A.neg;
   for(a = 0; a < L; a++)   C.pv[a] = C.hv[a] = 0;
   for(a = 0; a <= A.pl; a++)   C.pv[a] += A.pv[a];
   for(a = 0; a <= B.pl; a++)   C.pv[a] += B.pv[a];
   for(a = 0; a <= A.hl; a++)   C.hv[a] += A.hv[a];
   for(a = 0; a <= B.hl; a++) C.hv[a] += B.hv[a];
   h = (A.hl < B.hl) ? A.hl : B.hl;/*min*/
   for(a = h; a >= 1; a--)
       C.hv[a-1] += C.hv[a]/10, C.hv[a] %= 10;
   C.pv[0] += C.hv[0]/10, C.hv[0] %= 10;
   p = (A.hl > B.hl) ? A.hl : B.hl;/*max*/
   for(a = 0; a <= p; a++)
       C.pv[a+1] += C.pv[a]/10, C.pv[a] %= 10;
   p = L-1, h = L-1;
   while(C.pv[p] == 0 && p >= 0) p--;
   while(C.hv[h] == 0 && h >= 0) h--;
   C.pl = p, C.hl = h;
   return C;
}
BigNumber minus(BigNumber A, BigNumber B) {
   if(A.neg == 1 && B.neg == 1) {
       B.neg = 0, A.neg = 0;return minus(B, A);
   }
   else if(A.neg == 1 && B.neg == 0) {
       A.neg = 0;
       BigNumber C = addit(A, B);
       C.neg = 1;
       return C;
   }
   else if(A.neg == 0 && B.neg == 1) {
       B.neg = 0;return addit(A, B);
   }
   if(compare(A, B) == 0) {
       A.neg = 0;
       BigNumber C = minus(B, A);
       C.neg = 1;
       return C;
   }
   int a, p, h;
   BigNumber C;
   for(a = 0; a < L; a++)   C.pv[a] = C.hv[a] = 0;
   for(a = 0; a <= A.pl; a++) C.pv[a] = A.pv[a];
   for(a = 0; a <= A.hl; a++) C.hv[a] = A.hv[a];
   for(a = 0; a <= B.pl; a++) C.pv[a] -= B.pv[a];
   for(a = 0; a <= B.hl; a++) C.hv[a] -= B.hv[a];
   for(a = L-1; a >= 1; a--)
       while(C.hv[a] < 0)
           C.hv[a-1]--, C.hv[a] += 10;
   while(C.hv[0] < 0)
       C.pv[0]--, C.hv[0] += 10;
   for(a = 0; a < L; a++)
       while(C.pv[a] < 0)
           C.pv[a+1]--, C.pv[a] += 10;
   p = L-1, h = L-1;
   while(C.pv[p] == 0 && p >= 0) p--;
   while(C.hv[h] == 0 && h >= 0) h--;
   C.pl = p, C.hl = h;
   return C;
}
BigNumber multi(BigNumber A, BigNumber B) {
   BigNumber C;
   int Aa[2*L], Bb[2*L], Cc[4*L+1] = {};
   int a, b, Al = 0, Bl = 0, p, h;
   for(a = 0; a < L; a++)   C.pv[a] = C.hv[a] = 0;
   for(a = A.hl; a >= 0; a--)   Aa[Al++] = A.hv[a];
   for(a = 0; a <= A.pl; a++)   Aa[Al++] = A.pv[a];
   for(a = B.hl; a >= 0; a--)   Bb[Bl++] = B.hv[a];
   for(a = 0; a <= B.pl; a++)   Bb[Bl++] = B.pv[a];
   for(a = 0; a < Al; a++) {
       for(b = 0; b < Bl; b++)
           Cc[a+b] += Aa[a]*Bb[b];
   }
   for(a = 0; a < 2*L; a++)
       Cc[a+1] += Cc[a]/10, Cc[a] %= 10;
   int cut = (A.hl+1)+(B.hl+1);
   for(a = cut, b = 0; a < 4*L && b < L; a++, b++)
       C.pv[b] = Cc[a];
   for(a = cut-1, b = 0; a >= 0 && b < 100; a--, b++)
       C.hv[b] = Cc[a];
   p = L-1, h = L-1;
   while(C.pv[p] == 0 && p >= 0) p--;
   while(C.hv[h] == 0 && h >= 0) h--;
   C.pl = p, C.hl = h;
   C.pl = p, C.hl = h;
   C.neg = A.neg^B.neg;
   return C;
}
int Check(int a[], int b[], int lb, int lc) {
    int i;
    for(i = lb-1; i >= 0; i--) {
       if(a[i+lc] < b[i]) return 0;
       if(a[i+lc] > b[i]) return 1;
   }
   return 1;
}
BigNumber divis(BigNumber A, BigNumber B) {
   BigNumber C;
   int Aa[3*L] = {}, Bb[3*L] = {}, Cc[4*L+1] = {};
   int a, b, Al = L, Bl = 0, p, h, i, j;
   for(a = 0; a < L; a++)   C.pv[a] = C.hv[a] = 0;
   for(a = A.hl; a >= 0; a--)   Aa[Al++] = A.hv[a];
   for(a = 0; a <= A.pl; a++)   Aa[Al++] = A.pv[a];
   for(a = B.hl; a >= 0; a--)   Bb[Bl++] = B.hv[a];
   for(a = 0; a <= B.pl; a++)   Bb[Bl++] = B.pv[a];
   Al++, Bl++;
    for(i = Al-Bl; i >= 0; i--) {
       while(Check(Aa, Bb, Bl, i)){
           Cc[i]++;
            for(j = 0; j < Bl; j++) {
               Aa[i+j] -= Bb[j];
               if(Aa[i+j] < 0) {
                   Aa[i+j] += 10;
                   Aa[i+j+1]--;
                 }
           }
       }
    }
    int cut = L - B.hl + A.hl;
   for(a = cut, b = 0; a < 4*L && b < L; a++, b++)
       C.pv[b] = Cc[a];
   for(a = cut-1, b = 0; a >= 0 && b < L; a--, b++)
       C.hv[b] = Cc[a];
   p = L-1, h = L-1;
   while(C.pv[p] == 0 && p >= 0) p--;
   while(C.hv[h] == 0 && h >= 0) h--;
   C.pl = p, C.hl = h;
   C.neg = A.neg^B.neg;
   return C;
}
BigNumber msqrt(BigNumber A) {
   BigNumber C;
   int Aa[5*L] = {}, Cc[5*L+1] = {};
   int a, b, c, d, Al = 200 - A.hl - 1, Cl = 0, p, h;
   int cut;
   for(a = 0; a < L; a++)   C.pv[a] = C.hv[a] = 0;
   for(a = A.hl; a >= 0; a--)   Aa[Al++] = A.hv[a];
   for(a = 0; a <= A.pl; a++)   Aa[Al++] = A.pv[a];
    int B[5*L] = {}, top = Al, Bl = 1;
   for(a=(Al-1)/2*2 ; a >= 0; a -= 2, Bl++) {
        if(top > a+Bl-2)
           for(b = 9;b >= 0;b--) {
               int T[5*L]={0};
                B[0]=b;
               for(d=0;d <= Bl;d++) {
                    T[a+d] += (b*B[d]);
                    T[a+d+1] += (T[a+d]/10);
                    T[a+d] %= 10;
                }
               int find = 0;
                for(c = Bl+a+1; c >= 0; c--)
                   if(Aa[c] > T[c]) {find = 1;break;}
                    else if(T[c] > Aa[c]) break;
                if(find == 1 || c == -1) {
                   for(c = 0; c <= top;c++) {
                       Aa[c] -= T[c];
                       if(Aa[c]<0) {
                           int temp=((-Aa[c])/10+((-Aa[c])%10!=0));
                           Aa[c] += 10*temp;
                           Aa[c+1] -= temp;
                       }
                   }
                   break;
               }
            }
        else b = 0;
        Cc[Cl++] = b;
        if(a == 200)   cut = Cl;
        if(top != -1) {
           B[0] += b;
            for(b = 0; b <= Bl+1; b++)
               if(B[b] >= 10) {
                   B[b+1] += (B[b]/10);
                    B[b] %= 10;
                }
                else break;
           for(b = Bl+1;b >= 0; b--)   B[b+1] = B[b];
               B[0]=0;
            for(b = top;b >= 0; b--)
                if(Aa[b] != 0) {top = b;break;}
            if(b == -1) top = -1;
        }
    }
   for(a = cut-1, b = 0; a >= 0 && b < L; a--, b++)
       C.pv[b] = Cc[a];
   for(a = cut, b = 0; a < 3*L && b < L; a++, b++)
       C.hv[b] = Cc[a];
   p = L-1, h = L-1;
   while(C.pv[p] == 0 && p >= 0) p--;
   while(C.hv[h] == 0 && h >= 0) h--;
   C.pl = p, C.hl = h;
   C.neg = A.neg;
   return C;
}
void Print(BigNumber A) {
   int a;
   if(A.neg == 1)   printf("-");
   for(a = A.pl; a >= 0; a--)
       printf("%d", A.pv[a]);
   if(A.pl == -1) printf("0");
   if(A.hl != -1) printf(".");
   for(a = 0; a <= A.hl && a < 200; a++)
       printf("%d", A.hv[a]);

}
BigNumber Change(char s[]) {
    int a, l = strlen(s), cut = l, b, p, h, neg = 0;
    BigNumber C;
    if(s[0] == '-')    neg = 1;
    for(a = 0; a < L; a++)    C.pv[a] = C.hv[a] = 0;
    for(a = neg; a < l; a++)
        if(s[a] == '.') cut = a;
    for(a = cut-1, b = 0; a >= neg; a--, b++)
        C.pv[b] = s[a] - '0';
    for(a = cut+1, b = 0; a < l; a++, b++)
        C.hv[b] = s[a] - '0';
    p = L-1, h = L-1;
    while(C.pv[p] == 0 && p >= 0) p--;
    while(C.hv[h] == 0 && h >= 0) h--;
    C.pl = p, C.hl = h;
    C.neg = neg;
    return C;
}
main() {
    char x[200], y[20];
    while(scanf("%s", x) == 1) {
        BigNumber s1, s2, Ans, A;
        sprintf(y, "%d", 1);
        s1 = Change(x), s2 = Change(y);
        A = s1;
        sprintf(y, "%d", 0);
        Ans = Change(y);
        int a, t = 0;
        for(a = 1; a < 500; a += 2, t = 1-t) {
            if(t == 0)
                Ans = addit(Ans, divis(s1, s2));
            else
                Ans = minus(Ans, divis(s1, s2));
            /*Print(s1);
            puts("");
            Print(s2);
            puts("");*/
            s1 = multi(A, multi(s1, A));
            sprintf(y, "%d", (a+1)*(a+2));
            s2 = multi(s2, Change(y));

        }
        Print(Ans);
        puts("");
    }
    system("pause");
    return 0;
}
/*
-0.523593766666666666666666666666666666666666666666666666
*/

我要檢舉

#大數三角#sin

台長： Morris

您可能對以下文章有興趣

[機率][加分題77][Pascal分佈] Pascal random variable 證明

[中興資工][程設][作業HW13]

[其他題目] 由中根序列和后根序列重建二叉树

[練習] 打出垃圾碼，以 BrainFuck 為例

人氣(895) | 回應(0)| 推薦 (0)| 收藏 (0)| 轉寄
全站分類: 不分類 | 個人分類: 其他題目 |
此分類下一篇:目前世界最快－解數獨的程序, 個人實作
此分類上一篇:[手動][函式] 大數浮點數加減乘除開根號, 精準100位

回應(0)