ProjectEuler/c++/lib/lint.cpp
2013-04-17 14:34:39 +08:00

503 lines
10 KiB
C++
Raw Blame History

This file contains ambiguous Unicode characters

This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.

#include "lint.hpp"
void get(vector< complex<long double> > x)
{
for(uint i = x.size() - 1; i >= 0; i--) {
if(i > x.size()) break;
cout << x[i] << " ";
}
cout << endl << endl;
}
void get(vuint x)
{
for(uint i = x.size() - 1; i >= 0; i--) {
if(i > x.size()) break;
cout << x[i] << " ";
}
cout << endl << endl;
}
lint::lint() { }
lint::lint(string lintin)
{
operator =(lintin);
}
lint::~lint()
{ num.clear(); }
void lint::operator =(string strin)
{
sign = 1;
strin = del(strin);
if(strin[0] == '-') {
sign = -1;
strin = strin.substr(1, strin.size() - 1);
}
/* 整理输入使其全是数字并且首位无0 */
for(uint i = 0; i < strin.size(); i++)
if(strin[i] > '9' || strin[i] < '0') strin[i] = '0';
length = strin.size();
len = (length + unit - 1) / unit;
uint mod = length % unit;
if(num.capacity()) num.clear();
num.resize(len);
for(uint i = 0; i < len - 1; i++)
num[i] = atoi(strin.substr(length - (i + 1) * unit, unit).c_str());
num[len - 1] = atoi(strin.substr(0, (mod == 0) ? unit : mod).c_str());
}
/* 提取十进制数值 */
string lint::getd()
{ return vtos(num); }
vuint lint::getvec()
{
return num;
}
/* 提取 vector 中相应位置的数字 */
uint lint::operator [](uint locale)
{
if(locale >= len) return 0;
else return num[locale];
}
/* 提取 vector 的长度 */
uint lint::getlen()
{ return len; }
/* 提取整数长度 */
uint lint::getlength()
{ return length; }
int lint::getsign()
{ return sign; }
/* 定义加法 */
string lint::operator +(lint opnum)
{
if(sign * opnum.getsign() == -1) {
if(sign == -1) {
lint abs(getd());
return opnum - abs;
} else {
lint abs(opnum.getd());
lint absa(getd());
return absa - abs;
}
} else {
vuint sumint(1);
vuint sumshort(1);
if(len >= opnum.getlen()) {
sumint.resize(len);
sumint = getvec();
sumshort = opnum.getvec();
} else {
sumint.resize(opnum.getlen());
sumint = opnum.getvec();
sumshort = getvec();
}
sumint.resize(sumint.size() + 1);
sumint[sumint.size()] = 0;
uint temp;
for(uint i = 0; i < sumshort.size(); i++) {
temp = sumint[i] + sumshort[i];
sumint[i] = temp % cunit;
sumint[i + 1] += temp / cunit;
}
for(uint i = sumshort.size(); i < sumint.size() - 1; i++) {
temp = sumint[i];
sumint[i] = temp % cunit;
sumint[i + 1] += temp / cunit;
}
if(sign == -1) return "-" + vtos(sumint);
else return vtos(sumint);
}
}
string lint::operator -(lint opnum)
{
if(sign * opnum.getsign() == -1) {
if(sign == 1) return operator +(opnum);
else return "-" + operator +(opnum);
} else {
vuint difint(len >= opnum.getlen() ? len : opnum.getlen());
vuint small(len >= opnum.getlen() ? opnum.getlen() : len);
if(compare(opnum) >= 0) {
difint = num;
small = opnum.getvec();
} else {
difint = opnum.getvec();
small = num;
}
for(uint i = small.size() - 1; i >= 0; i--) {
if(i > small.size()) break;
if(difint[i] < small[i]) {
uint temp = i + 1;
while(difint[temp] == 0) {
difint[temp] = cunit - 1;
temp++;
if(temp >= difint.size()) break;
}
if(temp < difint.size()) difint[temp]--;
difint[i] += cunit;
}
difint[i] -= small[i];
}
if(sign * compare(opnum) == 1) return vtos(difint);
else return "-" + vtos(difint);
}
}
/* 定义乘法,用 FFT 实现 */
string lint::operator *(lint opnum)
{
uint a_size = len, b_size = opnum.getlen();
vector< complex<long double> > a(a_size), b(b_size);
for(uint i = 0; i < a_size; i++) a[i].real() = num[i];
for(uint i = 0; i < b_size; i++) b[i].real() = (opnum[i]);
uint n;
a_size += b_size;
for (n = a_size; n != (n&-n); n += (n&-n));
a.resize(n);
b.resize(n);
FFT(a, n);
FFT(b, n);
for (uint i = 0; i < n; i++)
a[i] *= b[i];
IFFT(a, n);
for (uint i = 0; i < a_size-1; i++) {
a[i + 1] += a[i].real() / 1000.;
a[i] = fmodl(a[i].real(), 1000.L);
a[i].real() += eps;
}
vuint mulout(n);
for(uint i = 0; i < a.size(); i++) mulout[i] = a[i].real();
if(sign * opnum.getsign() == -1) return "-" + vtos(mulout);
else return vtos(mulout);
}
int lint::operator ==(lint opnum)
{
if(sign * opnum.getsign() > 0) {
if(sign == 1) return compare(opnum);
if(sign == -1) return -(compare(opnum));
else return 0;
} else {
if(sign == 1) return 1;
else return -1;
}
}
int lint::compare(lint opnum)
{
if(length > opnum.getlength()) return 1;
if(length < opnum.getlength()) return -1;
else {
for(uint i = len - 1; i >= 0; i--) {
if(i > len - 1) break;
if(num[i] > (opnum[i])) return 1;
if(num[i] < (opnum[i])) return -1;
}
return 0;
}
}
/* 定义除法 */
string lint::operator /(lint opnum)
{
if(compare(opnum) == -1) return "0";
if(compare(opnum) == 0) return "1";
else {
uint move = len - opnum.getlen();
vuint devint(move + 1);
vuint open = getvec();
open.resize(len + 1);
open[len] = 0;
for(uint i = move; i >= 0; i--) {
if(i > move) break;
uint high = i + opnum.getlen() - 1;
if(open.capacity() <= high) {
open.resize(high + 2);
open[high + 1] = 0;
}
uint temp = (open[high] + open[high + 1] * cunit) / opnum.getvec().back() + 2;
vuint tryvec;
int count = 1;
do {
if(--temp > cunit * cunit) temp = 1;
tryvec = vmul(opnum.getvec(), temp, i);
count++;
} while(vless(open, tryvec) && count < 20);
devint[i] = temp;
open = vdif(open, tryvec);
open[open.size()] = 0;
tryvec.clear();
}
if(sign * opnum.getsign() == 1) return vtos(devint);
else return "-" + vtos(devint);
}
}
string lint::operator %(lint opnum)
{
if(sign == 1 && compare(opnum) == -1) return getd();
else {
lint div = getd();
lint q = (div / opnum);
lint mod = (div - (q * opnum));
if(sign == 1) return mod.getd();
else return (opnum - mod);
}
}
vuint lint::getbit()
{
vuint out(len * 3);
lint opnum = getd();
lint bit("2");
uint i;
for(i = 0; i < len * 3; i++) {
if(opnum % bit == "1") out[i] = 1;
else out[i] = 0;
opnum = opnum / bit;
}
while(opnum.getd() != "0") {
out.resize(i + 2);
if(opnum % bit == "1") out[i] = 1;
else out[i] = 0;
opnum = opnum / bit;
i++;
}
return del(out);
}
string lint::operator ^(lint opnum)
{
lint power("1");
lint base = getd();
for(uint i = opnum.getbit().size() - 1; i >= 0; i--) {
if(i > opnum.getbit().size()) break;
power = (power * power);
if(opnum.getbit()[i] == 1) power = (power * base);
}
return del(power.getd());
}
/* 消去首位0 */
string del(string undel)
{
uint flag = 0;
while(undel.substr(0, 2) == "--") undel = undel.substr(2, undel.size() - 2);
if(undel[0] == '-') {
undel = undel.substr(1, undel.size() - 1);
flag = 1;
}
while(undel.size() > 1 && undel[0] == '0')
undel = undel.substr(1, undel.size() - 1);
if(undel.size() == 1 && undel[0] == '0') return "0";
if(flag) return "-" + undel;
return undel;
}
vuint del(vuint undel)
{
while(undel[undel.size() - 1] == 0) undel.pop_back();
vuint out(undel.size());
for(uint i = 0; i < undel.size(); i++) out[i] = undel[i];
return out;
}
/* 整数到字符串的转换 */
string itos(uint intin)
{
string strout = "";
if(intin == 0) return "0";
while(intin > 0) {
char temp = intin % carry + '0';
strout = temp + strout;
intin /= carry;
}
return strout;
}
/* 向量到字符串的转换 */
string vtos(vuint vin)
{
string strout = "";
string temp = "";
for(uint i = vin.size() - 1; i >= 0; i--) {
if(i > vin.size() - 1) break;
temp = itos(vin[i] % cunit);
while(temp.size() < unit) temp = "0" + temp;
strout += temp;
}
return del(strout);
}
/**
向量转比特流
string vtobit(vuint vin)
{
string strout = "";
for(uint i = vin.size() - 1; i >= 0; i--) {
if(i > vin.size()) break;
if(vin[i]) strout += "1";
else strout += "0";
}
return del(strout);
}*/
/* 向量和整数的乘法 */
vuint vmul(vuint op, uint opnum, uint power = 0)
{
if(power) {
op.resize(op.size() + power);
for(uint i = op.size() - 1; i >= power; i--) {
if(i > op.size()) break;
op[i] = op[i - power];
}
for(uint i = power - 1; i >= 0; i--) {
if(i > op.size()) break;
op[i] = 0;
}
}
lint big(vtos(op));
lint small(itos(opnum));
lint out(big * small);
return out.getvec();
}
/* 向量减法 */
vuint vdif(vuint big, vuint small)
{
lint bigint(vtos(big));
lint smallint(vtos(small));
if(bigint.compare(smallint) == -1) return big;
else {
lint out(bigint - smallint);
return out.getvec();
}
}
/* 比较大小 */
uint vless(vuint big, vuint small)
{
lint bigint(vtos(big));
lint smallint(vtos(small));
if(bigint.compare(smallint) == -1) return 1;
else return 0;
}
/* 快速傅立叶变换 */
void FFT(vector< complex<long double> > &x, uint n)
{
uint i,j,k,t;
for (i = 0; i < n; ++i) {
j = 0;
for (t = i, k = n; k /= 2; t /= 2)
j = j*2+t%2;
if (j > i) swap(x[j], x[i]);
}
for (k = 2; k <= n; k *= 2) {
const complex<long double> omega_unit(cosl(2*PI/k), sinl(2*PI/k));
for (i = 0; i < n; i += k) {
complex<long double> omega(1, 0);
for (j = 0; j < k/2; ++j) {
complex<long double> t = omega*x[i+j+k/2];
x[i+j+k/2] = x[i+j]-t;
x[i+j] += t;
omega *= omega_unit;
}
}
}
}
void IFFT(vector< complex<long double> > &x, uint n)
{ uint i,j,k,t;
for (i = 0; i < n; ++i) {
j = 0;
for (t = i, k = n; k /= 2; t /= 2)
j = j*2+t%2;
if (j > i) swap(x[j], x[i]);
}
for (k = 2; k <= n; k *= 2) {
const complex<long double> omega_unit(cosl(-2*PI/k), sinl(-2*PI/k));
for (i = 0; i < n; i += k) {
complex<long double> omega(1, 0);
for (j = 0; j < k/2; ++j) {
complex<long double> t = omega*x[i+j+k/2];
x[i+j+k/2] = x[i+j]-t;
x[i+j] += t;
omega *= omega_unit;
}
}
}
for (i = 0; i < n; ++i)
x[i] /= n;
}
/* 求两数最大公约数 */
string gcd(lint a, lint b)
{
if(b.getd() == "0") return a.getd();
else {
lint newb = a % b;
return gcd(b, newb);
}
}