503 lines
10 KiB
C++
503 lines
10 KiB
C++
#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);
|
||
}
|
||
}
|