|
Реализация вещественной длинной арифметики (с использованием класса vector)
Необходимо реализовать операции сложения, вычитания, умножения, деления и нахождения корня для чисел большой длины с плавающей точкой.
Листинг C++
Листинг C++
#define forp(i, a, b, c) for (i = (a); (c) > 0 ? i <= (b) : i >= (b); i += (c))
#define forn(i, n) forp (i, 0, (n) - 1, 1)
#define form(i, n) forp (i, (n) - 1, 0, - 1)
#define forv(it, v) for (it = (v).begin (); it != (v).end (); ++ it)
#define uint long long
//--------------------- Defines ------------------------------------------------
#define sint int
#define vsint vector < sint >
#define sz size ()
#define at(v, i) (0 <= (i) && (i) < (v).size () ? v[i] : 0)
#define erase_null(v) while ((v).size () > 1 && (v).back () == 0) (v).pop_back ();
#define next_mod(v, i) if (0 <= (i) + 1 && (i) + 1 < (v).size ()) { v[(i) + 1] += v[(i)] / 10; v[(i)] %= 10; }
//--------------------- Shifts, Compare ----------------------------------------
vsint shift_left (vsint s, int k)
{
if (k <= 0) return s;
s.insert (s.begin (), k, 0);
return s;
}
vsint shift_right (vsint s, int k)
{
if (k <= 0) return s;
s.erase (s.begin (), s.begin () + min (k, (int)s.sz));
return s;
}
int cmp (vsint &a, vsint &b)
{
for (int i = max (a.sz, b.sz) + 1; i >= 0; -- i)
{
if (at (a, i) > at (b, i)) return 1;
if (at (a, i) < at (b, i)) return - 1;
}
return 0;
}
//--------------------- Add, Sub, Mul, Div -------------------------------------
vsint add (vsint &a, vsint &b)
{
vsint c (max (a.sz, b.sz) + 1, 0);
for (int i = 0; i < c.sz - 1; ++ i)
{
c[i] += at (a, i) + at (b, i);
next_mod (c, i);
}
erase_null (c);
return c;
}
vsint sub (vsint &a, vsint &b)
{
vsint c (a.sz, 0);
for (int i = 0; i < c.sz; ++ i)
{
c[i] += at (a, i) - at (b, i);
if (c[i] < 0 && i + 1 < c.sz)
{
c[i] += 10;
-- c[i + 1];
}
}
erase_null (c);
return c;
}
void sub_from (vsint &a, vsint &b)
{
for (int i = 0; i < a.sz; ++ i)
{
a[i] -= at (b, i);
if (a[i] < 0 && i + 1 < a.sz)
{
a[i] += 10;
-- a[i + 1];
}
}
erase_null (a);
}
vsint mul (vsint &a, vsint &b)
{
vsint c (a.sz + b.sz + 2, 0);
for (int i = 0; i < a.sz; ++ i)
for (int j = 0; j < b.sz; ++ j)
{
c[i + j] += a[i] * b[j];
next_mod (c, i + j);
}
for (int i = 0; i < c.sz; ++ i)
next_mod (c, i);
erase_null (c);
return c;
}
vsint div (vsint &a, vsint &b, vsint &r)
{
r.assign (1, 0);
vsint s (a.sz, 0);
for (int i = a.sz - 1; i >= 0; -- i)
{
r.insert (r.begin (), a[i]);
while (cmp (r, b) >= 0)
{
++ s[i];
sub_from (r, b);
}
}
erase_null (s);
erase_null (r);
return s;
}
//--------------------- Algo ---------------------------------------------------
vsint sqrt (vsint &a)
{
int n = (a.size () + 1) / 2;
vsint x (n, 0);
for (int i = n - 1; i >= 0; -- i)
for (x[i] = 9; x[i] > 0; -- x[i])
{
vsint xx = mul (x, x);
if (cmp (xx, a) <= 0)
break;
}
return x;
}
//------------------------------------------------------------------------------
//--------------------- BigDouble ----------------------------------------------
const int DoubleSize = 32;
class bdouble
{
public :
vsint a;
char sign;
bdouble (double x);
void sscan (string s);
void normal ();
void print ();
friend bdouble operator - (bdouble a);
friend int cmp (bdouble a, bdouble b);
friend bdouble operator + (bdouble a, bdouble b);
friend bdouble operator - (bdouble a, bdouble b);
friend bdouble operator * (bdouble a, bdouble b);
friend bdouble operator / (bdouble a, bdouble b);
};
bdouble :: bdouble (double x = 0.0)
{
char s[100] = {0};
sprintf (s, "%.15lf", x);
sscan (s);
}
void bdouble :: sscan (string s)
{
sign = 1;
if (s[0] == '-')
{
sign = - 1;
s = s.substr (1);
}
int i, j;
int n = s.size ();
int z = s.find ('.');
if (z == string :: npos) z = n;
a.assign (DoubleSize, 0);
for (i = z + 1; i < n; ++ i)
a[DoubleSize - (i - (z + 1)) - 1] = s[i] - '0';
for (i = z - 1; i >= 0; -- i)
a.push_back (s[i] - '0');
}
void bdouble :: normal ()
{
while (a.sz < DoubleSize + 1)
a.push_back (0);
}
void bdouble :: print ()
{
if (sign < 0)
printf ("-");
for (int i = a.sz - 1; i >= 0; -- i)
{
printf ("%d", a[i]);
if (i == DoubleSize)
printf (".");
}
}
bdouble operator - (bdouble a)
{
a.sign = - a.sign;
return a;
}
int cmp (bdouble a, bdouble b)
{
if (a.sign > b.sign) return 1;
if (a.sign < b.sign) return - 1;
return cmp (a.a, b.a) * a.sign;
}
bdouble operator + (bdouble a, bdouble b)
{
if (a.sign != b.sign)
if (b.sign < 0) return a - (- b);
else return b - (- a);
bdouble c;
c.a = add (a.a, b.a);
c.sign = a.sign;
c.normal ();
return c;
}
bdouble operator - (bdouble a, bdouble b)
{
if (cmp (a, b) < 0) return - (b - a);
if (b.sign < 0) return a + (- b);
bdouble c;
c.a = sub (a.a, b.a);
c.sign = 1;
c.normal ();
return c;
}
bdouble operator * (bdouble a, bdouble b)
{
bdouble c;
c.a = mul (a.a, b.a);
c.a = shift_right (c.a, DoubleSize);
c.sign = a.sign * b.sign;
c.normal ();
return c;
}
bdouble operator / (bdouble a, bdouble b)
{
bdouble c;
a.a = shift_left (a.a, DoubleSize);
vsint r;
c.a = div (a.a, b.a, r);
c.sign = a.sign * b.sign;
c.normal ();
return c;
}
//--------------------- Algo ---------------------------------------------------
bdouble bsqrt (bdouble a)
{
bdouble x;
a.a = shift_left (a.a, DoubleSize);
x.a = sqrt (a.a);
x.sign = 1;
x.normal ();
return x;
}
15:44
06.08.2009
По всем вопросам обращаться: rumterg@gmail.com