郭先強演算法

❶ 徵求n階乘的優化演算法

這是被我初學時寫的如今已我拋棄的"演算法", 其實也不是什麼演算法, 只不過類設計一個, 而且類的實現還沒有完善, 因為寫到一半發現我的構思有多麼垃圾我就沒再寫下去了, 支持4294967295位有符號和無符號整數, 你可以把UBI_SIGNED給undef了對於無符號整型可能要快一點, 我也沒算過100000!有沒有超過這個范圍, 你要算的話自己去算下看看, 我沒那個時間.

UBI意為unbounded integer

忘了說你可以去這個網站看下一個專門為理論上的"無限大"數做的的C++庫:

http://www.apfloat.org/apfloat/

apfloat的意思是arbitrary precision float, 這個人寫的演算法很吊, 只能這樣說, 2.6億位的π只用幾乎1秒多就出來了, 算100W!也不在話下, 用的大約是數論變換, 你要演算法的話可以去參考下, 另:
其實the art of computer programming vol2 siminumerical algorithm有很多關於這方面的介紹的, 我曾經瞥過幾眼, 但由於本人數學太爛, 看懂幾乎是 impossible -_-~~. 所以也只好作罷...

我寫這個爛代碼計算乘法很"慢", 但計算加法還馬虎, 我曾經用跌代法算Fibonacci的第10W項貌似也只用了很短的時間(感覺不長, 具體記不得了).

計算了一個1000!, 3秒左右:

1239862902
9445909974
1418278094
5573543251
2912073791
0878297308
8932076716
3650241536
0821333186
4516076535
6153071277
8114194545
3408243920
9880051954
5832036786
6848259012
3448259932
8558600301
2188525247
4228458623
5904339901
0164192106
0241866493
1399694290
986355
8633969099
2154399457
8184009699
0277534720
0000000000
0000000000
0000000000
00000000

代碼:

#include <iostream>
#include <vector>
#include <string>
#include <sstream>
#include <stdexcept>
using namespace std;

#define UBI_SIGNED

class UBI
{
typedef unsigned char byte;
typedef signed char diff_byte;
typedef vector<byte> mem_type;
typedef mem_type::size_type size_type;

#define SMALL(result) result.first
#define BIG(result) result.second
typedef pair<const UBI*, const UBI*> cmp_result;
static inline byte HiByte(byte b)
{
return b >> 4;
}
static inline byte LoByte(byte b)
{
return b & 0xF;
}
static inline byte MakeByte(byte a, byte b)
{
return (a << 4) | (b & 0xF);
}

void Init(string str);
short Strcmp(const string& lhs, const string& rhs)const;
cmp_result UBICmp(const UBI& lhs, const UBI& rhs)const;

#ifndef UBI_SIGNED
void CheckVal(const string& str)const
{
if(str[0] == '-')
throw domain_error
("Try to initial or assign an negtive value" \
"to a number that is positive.");
}
#endif

public:
UBI(size_type digit = 20);
UBI(const char* c_str, size_type digit = 20);
UBI(const string& str, size_type digit = 20);
UBI(const UBI& UBI);

UBI operator + (const UBI& rhs)const;
UBI operator - (const UBI& rhs)const;
UBI operator * (const UBI& rhs)const;
UBI operator / (const UBI& rhs)const;

const UBI& operator = (const UBI& rhs);
const UBI& operator += (const UBI& rhs);
const UBI& operator -= (const UBI& rhs);
const UBI& operator *= (const UBI& rhs);
const UBI& operator /= (const UBI& rhs);

string ToString()const;
void clear();

private:
mem_type mem;
#ifdef UBI_SIGNED
bool isSigned;
bool IsSigned()const;
mutable bool mutex;
#endif
};

UBI::UBI(size_type digit)
{
#ifdef UBI_SIGNED
isSigned = false;
mutex = false;
#endif
mem.reserve(digit);
}

void UBI::Init(string str)
{
if(str[0] == '-')
str.erase(0, 1);
int i;
for(i = str.size() - 1; i-1 >= 0; i -= 2)
{
mem.push_back(MakeByte(LoByte(str[i]), LoByte(str[i-1])));
}

if(i == 0)mem.push_back(MakeByte(LoByte(str[0]), 0));
}

short UBI::Strcmp(const string& lhs, const string& rhs)const
{
if(lhs.size() > rhs.size())
return 1;
else if(lhs.size() < rhs.size())
return -1;
else if(lhs > rhs)
return 1;
else if(lhs < rhs)
return -1;
return 0;
}

UBI::cmp_result UBI::UBICmp(const UBI& lhs, const UBI& rhs)const
{
string sLhs(lhs.ToString());
string sRhs(rhs.ToString());

if(Strcmp(sLhs, sRhs) == 0)
return make_pair(&lhs, &lhs);

#ifdef UBI_SIGNED
if(sLhs[0] == '-')sLhs.erase(0, 1);
if(sRhs[0] == '-')sRhs.erase(0, 1);
#endif
const UBI* small =
Strcmp(sLhs, sRhs) == -1 ?
this : &rhs;

const UBI* big =
&rhs == small ?
this : &rhs;

return make_pair(small, big);
}

#ifdef UBI_SIGNED
bool inline UBI::IsSigned()const
{
return isSigned;
}
#endif

UBI::UBI(const char* c_str, size_type digit)
{
#ifdef UBI_SIGNED
isSigned = *c_str == '-' ? true : false;
mutex = false;
#else
CheckVal(string(c_str));
#endif
if(string(c_str).size() > digit)
mem.reserve(string(c_str).size());
else
mem.reserve(digit);
Init(string(c_str));
}

UBI::UBI(const string& str, size_type digit)
{
#ifdef UBI_SIGNED
isSigned = str[0] == '-' ? true : false;
mutex = false;
#else
CheckVal(str);
#endif
if(str.size() > digit)
mem.reserve(str.size());
else
mem.reserve(digit);
Init(str);
}

UBI::UBI(const UBI& UBI) : mem(UBI.mem)
{
#ifdef UBI_SIGNED
isSigned = UBI.isSigned;
mutex = false;
#endif
}

const UBI& UBI::operator = (const UBI& rhs)
{
if(this == &rhs)return *this;

mem.assign(rhs.mem.begin(), rhs.mem.end());
#ifdef UBN_SIGNED
isSigned = rhs.isSigned;
#endif
return *this;
}

UBI UBI::operator + (const UBI& rhs)const
{
if(mem.empty() || mem.size() == 1 && mem[0] == 0x0)
return rhs;
if(rhs.mem.empty() || rhs.mem.size() == 1 && rhs.mem[0] == 0x0)
return *this;

#ifdef UBI_SIGNED
if(!mutex)
{
if(isSigned && !rhs.isSigned || !isSigned && rhs.isSigned)
{
mutex = true;
return *this - rhs;
}
}else
mutex = false;
#endif

cmp_result result(UBICmp(*this, rhs));

byte prevHiRemain = 0;
size_type smallSize = SMALL(result)->mem.size();
UBI tempUBI;
for(size_type i = 0; i < BIG(result)->mem.size(); ++i)
{
byte tempHi =
HiByte(BIG(result)->mem[i]) +
(i < smallSize ? HiByte(SMALL(result)->mem[i]) : 0)
+ prevHiRemain;

byte tempLo =
LoByte(BIG(result)->mem[i]) +
(i < smallSize ? LoByte(SMALL(result)->mem[i]) : 0);

prevHiRemain = 0;
if(tempHi > 9)
{
tempLo += tempHi / 10;
tempHi = tempHi % 10;
}
if(tempLo > 9)
{
prevHiRemain = tempLo / 10;
tempLo = tempLo % 10;
}
tempUBI.mem.push_back(MakeByte(tempHi, tempLo));
}
if(prevHiRemain != 0)
tempUBI.mem.push_back(MakeByte(prevHiRemain, 0));

#ifdef UBI_SIGNED
tempUBI.isSigned = this->isSigned ? true : false;
#endif
return tempUBI;
}

const UBI& UBI::operator += (const UBI& rhs)
{
return *this = *this + rhs;
}

UBI UBI::operator - (const UBI& rhs)const
{
#ifdef UBI_SIGNED
if(!mutex)
{
if(!isSigned && rhs.isSigned || isSigned && !rhs.isSigned)
{
mutex = true;
return *this + rhs;
}
}else
mutex = false;

#endif

cmp_result result(UBICmp(*this, rhs));
if(SMALL(result) == BIG(result))
return UBI("0");

byte prevHiBorrow = 0;
size_type smallSize = SMALL(result)->mem.size();
UBI tempUBI;
for(size_type i = 0; i < BIG(result)->mem.size(); ++i)
{
diff_byte tempHi =
HiByte(BIG(result)->mem[i]) -
(i < smallSize ? HiByte(SMALL(result)->mem[i]) : 0)
- prevHiBorrow;
diff_byte tempLo =
LoByte(BIG(result)->mem[i]) -
(i < smallSize ? LoByte(SMALL(result)->mem[i]) : 0);

prevHiBorrow = 0;
if(tempHi < 0)
{
tempLo -= 1;
tempHi = 10 + tempHi;
}
if(tempLo < 0)
{
prevHiBorrow = 1;
tempLo += 10;
}
tempUBI.mem.push_back(MakeByte(tempHi, tempLo));
}

#ifdef UBI_SIGNED
if(this == BIG(result) && this->isSigned ||
this == SMALL(result) && !this->isSigned)
tempUBI.isSigned = true;
#endif

return tempUBI;
}

UBI UBI::operator * (const UBI& rhs)const
{
if(this->mem[mem.size()-1] == 0 || rhs.mem[rhs.mem.size()-1] == 0)
return UBI("0");
if(this->mem.size() == 1 && this->mem[0] == 0x10)
return rhs;
if(rhs.mem.size() == 1 && rhs.mem[0] == 0x10)
return *this;

cmp_result result(UBICmp(*this, rhs));

UBI tempUBI;
for(size_type i = 0, k = 0; i < SMALL(result)->mem.size(); ++i)
{
byte SmHi = HiByte(SMALL(result)->mem[i]);
byte SmLo = LoByte(SMALL(result)->mem[i]);
byte prevLoRemain = 0;

UBI Hi;
for(size_type j = 0; j < BIG(result)->mem.size(); ++j)
{
byte tempHi = SmHi * HiByte(BIG(result)->mem[j]) + prevLoRemain;
byte tempLo = SmHi * LoByte(BIG(result)->mem[j]);
prevLoRemain = 0;

if(tempHi > 9)
{
tempLo += tempHi / 10;
tempHi %= 10;
}
if(tempLo > 9)
{
prevLoRemain = tempLo / 10;
tempLo %= 10;
}
Hi.mem.push_back(MakeByte(tempHi, tempLo));
}
if(prevLoRemain != 0)
{
Hi.mem.push_back(MakeByte(prevLoRemain, 0));
prevLoRemain = 0;
}
if(k > 0)
{
string _10(Hi.ToString());
_10.resize(_10.size()+k, '0');
Hi = UBI(_10);
}
++k;
tempUBI += Hi;

UBI Lo;
for(size_type j = 0; j < BIG(result)->mem.size(); ++j)
{
byte tempHi = SmLo * HiByte(BIG(result)->mem[j]) + prevLoRemain;
byte tempLo = SmLo * LoByte(BIG(result)->mem[j]);
prevLoRemain = 0;

if(tempHi > 9)
{
tempLo += tempHi / 10;
tempHi %= 10;
}
if(tempLo > 9)
{
prevLoRemain = tempLo / 10;
tempLo %= 10;
}
Lo.mem.push_back(MakeByte(tempHi, tempLo));
}
if(prevLoRemain != 0)
{
Lo.mem.push_back(MakeByte(prevLoRemain, 0));
prevLoRemain = 0;
}
if(k > 0)
{
string _10(Lo.ToString());
_10.resize(_10.size()+k, '0');
Lo = UBI(_10);
}
++k;
tempUBI += Lo;
}
return tempUBI;
}

string UBI::ToString()const
{
if(mem.empty() || mem.size() == 1 && mem[0] == 0)
return string("0");
string temp;
temp.reserve(mem.size()*2);

#ifdef UBI_SIGNED
if(isSigned)
temp.push_back('-');
#endif
int i = mem.size()-1;
while(mem[i] == 0 && i > 0)--i;

for(; i >= 0; --i)
{
temp.push_back(LoByte(mem[i]) + '0');
temp.push_back(HiByte(mem[i]) + '0');
}

unsigned zero = temp.find_first_of("0");
#ifdef UBI_SIGNED
if(zero == 0 || temp[0] == '-' && zero == 1)
temp.erase(zero, 1);
#else
if(zero == 0)
temp.erase(zero, 1);
#endif
return temp;
}

void UBI::clear()
{
#ifdef UBI_SIGNED
isSigned = false;
#endif
mem.clear();
}

template <class T, class U>
T lexical_cast(U u)
{
stringstream sstrm;
sstrm << u;
T t;
sstrm >> t;
return t;
}

int main()
{
UBI a("1");
for(int i = 2; i <= 1000; ++i)
{
a = a * lexical_cast<string>(i);
}
cout << a.ToString();
}