Длинная вещественная арифметика. Пишем свои вещественные числа на C++

У вас бывало так, что при попытке подсчитать какую-нибудь математическую константу с достаточно большой точностью (например, с помощью вычисления суммы ряда), вам этого всё равно не удавалось? Это происходит из-за того, что встроенные типы данных для выполнения операций с вещественными числами (float, double, long double) имеют свои ограничения как на минимальные с максимальными значения, так и на количество допустимых знаков после запятой.

Так что, с их помощью не получится вычислить число e с точностью до 100 знака? Однозначно нет. Используя встроенные типы чисел с плавающей точкой, экспоненту до сотого десятичного знака вычислить не удастся. А вот с помощью других вещественных чисел, в которых бы, например, отсутствовали (в разумных пределах конечно же) всякие ограничения на возможный диапазон значений и знаков после запятой эта задача легко решилась!

Вот мы и подошли к основной задаче: требуется создать такой тип данных, который может хранить вещественные числа неограниченной длины и выполнять с ними базовые арифметические операции (+, -, *, /). В качестве языка программирования, как легко догадаться из названия статьи, выберем C++.

Для самых нетерпеливых и жаждущих код, добро пожаловать на github проекта.

Как хранить вещественные числа неограниченной длины?

Будем хранить число в виде следующей структуры:

• вектор цифр (digits): игнорируем запятую, просто записывая все цифры числа
• знак (sign): число -1 или 1 для соответственно отрицательных и положительных чисел
• экспонента (exponent): число, на которое нужно свдинуть цифры так, чтобы запятая оказалась в исходном месте.

Таким образом, мы будем представлять числа в виде sign * 0.d1d2...dn * 10exponent.

Если, например, мы захотим сохранить число 1234.5678 в нашу структуру, то это будет выглядеть следующим образом: sign: 1, digits: 12345678, exponent: 4

А если потребуется сохранить число -0.00012345678, то получим такую структуру: sign: -1, digits: 12345678, exponent: -3

Примечания:

• Эспонента фактически показывает, сколько цифр находится в целой части (или сколько находится нулей после запятой, если экспонента отрицательная)
• Поскольку цифры будут храниться в векторе, то длина чисел всё-таки будет ограничена числом 2⁶⁴ (максимальное количество элементов в векторе), однако это число настолько велико, что вряд ли вам когда-либо удастся его достичь.

Создаём класс длинных вещественных чисел

Создадим класс с конструктором по-умолчанию, конструктором из строки (string) и вышеперечисленными закрытыми полями:

#include <iostream>
#include <string>
#include <vector>

using namespace std;

class LongDouble {
	int sign; // знак числа
	vector<int> digits; // вектор цифр
	long exponent; // экспонента для свдига

public:
	LongDouble(); // конструктор по-умолчанию
	LongDouble(const LongDouble& x); // конструктор копирования
	LongDouble(const string &value); // конструктор из строки
};

В конструкторе по-умолчанию будем создавать нулевое число, то есть в векторе будет всего одна цифра - 0, а экспонента будет равна 1, так как 0.000 имеет ровно одну цифру в целой части.

// конструктор по-умолчанию
LongDouble::LongDouble() {
	sign = 1; // ноль будем считать положительным числом
	digits = vector<int>(1, 0); // создаём вектор из 1 нулевой цифры
	exponent = 1;
}

// конструктор копирования
LongDouble::LongDouble(const LongDouble& x) {
	sign = x.sign; // копируем знак
	exponent = x.exponent; // копируем экспоненту
	digits = vector<int>(x.digits); // копируем цифры
}

// конструктор из строки
LongDouble::LongDouble(const string& value)
	size_t index;

	// если первый символ строки - минус, значит число отрицательное
	if (value[0] == '-') {
		sign = -1; // знак отрицательного числа -1
		index = 1; // начинать идти нужно будет с первого символа
	} else {
		sign = 1; // иначе число положительное
		index = 0; // и идти нужно будет с нулевого символа 
	}

	exponent = value.length() - index; // предполагаем, что всё число будет целым, а значит и экспонента будет равна длине строки

	// идём по всей строке
	while (index < value.length()) {
		if (value[index] == '.')  // если встретили точку
			exponent = sign == 1 ? index : index - 1; // запоминаем экспоненту, так как целая часть закончилась
		else
			digits.push_back(value[index] - '0'); // иначе вставляем в вектор очередную цифру

		index++; // переходим к новому символу
	}
}

Стоит отметить, что только ноль возможно записать различными способами, так как при разных экспонентах, он всё равно останется нулём.

Может оказаться такая ситуация, что в строке будут лишние нули как с правой, так и с левой сторон. Для этого создадим в классе приватную функцию void removeZeroes(), которая уберёт лишние нули справа и слева. Причём убирая нули слева (после запятой), нужно учитывать, что это приводит к увеличению числа в 10 раз, а значит нужно уменьшать экспоненту на 1 при каждом таком удалении.

void LongDouble::removeZeroes() {
	size_t n = max((long) 1, exponent); // определяем, до какого момента нужно будет идти для удаления нулей справа

	// пока справа нули
	while (digits.size() > n && digits[digits.size() - 1] == 0)
		digits.erase(digits.end() - 1); // удаляем их

	// пока число цифр больше одной и вектор начинается с нуля
	while (digits.size() > 1 && digits[0] == 0) {
		digits.erase(digits.begin()); // удаляем нули из начала вектора
		exponent--; // и уменьшаем экспоненту
	}

	// если справа всё ещё остались нули
	while (digits.size() > 1 && digits[digits.size() - 1] == 0)
		digits.erase(digits.end() - 1); // то удалим их

	// если в результате осталась всего одна нулевая цифра, то превратим её в честный ноль
	if (digits.size() == 1 && digits[0] == 0) {
		exponent = 1;
		sign = 1;
	}
}

Чтобы корректно хранить число, полученное из строки, будем вызывать написанную функцию в самом конце конструктора из строки

Добавляем вывод объектов

Теперь наш класс способен хранить числа неограниченной длины, получая их из строки. Хочется проверить, что всё сделано правильно, но для этого потребуется написать вывод этого объекта. Выводить объект будем в стандарнтый поток вывода ostream с помощью перегрузки оператора <<:

ostream& operator<<(ostream& os, const LongDouble& value) {
	if (value.sign == -1)
		os << '-'; // если число отрицательное, то сначала выведем знак минус

	// если экспонента положительна, то у числа ненулевая целая часть
	if (value.exponent > 0) {
		size_t i = 0;
		size_t e = value.exponent;

		// выводим первые exponent цифр (или все цифры, если экспонента больше) числа чтобы вывести целую часть
		while(i < value.digits.size() && i < e)
			os << value.digits[i++];

		// если экспонента больше цифр числа, то выводим нули, чтобы дойти до экспоненты
		while (i < e) {
			os << "0";
			i++;
		}

		// если цифры ещё остались
		if (i < value.digits.size()) {
			os << "."; // то выводим точку

			// и выводим оставшиеся цифры как дробную часть
			while(i < value.digits.size())
				os << value.digits[i++];
		}
	}
	else { // иначе эспонента отрицательна или нулевая
		os << "0."; // выводим нулевую целую часть с точкой

		// выводим |exponent| нулей (если экспонента нулевая, то не будет ни одного нуля)
		for (long i = 0; i < -value.exponent; i++)
			os << "0";

		// выводим все цифры числа
		for (size_t i = 0; i < value.digits.size(); i++)
			os << value.digits[i];
	}

	return os; // возвращаем поток вывода
}

Отлично! Теперь мы не только можем создать длинное вещественное число из строки, но и вывести его. А это уже весьма немало! Теперь никто не помешает создать число 123154654.1234543245643245643456434565434567543234567876543234567, которое при выводе останется абсолютно тем же самым, а не превратится в "обрубок" из-за округления встроенных вещественных типов языка C++.

Чтобы проверить вышесказанное, напишем основную программу, в которой и создадим объект этого числа из строки:

int main() {
	LongDouble a("123154654.1234543245643245643456434565434567543234567876543234567"); // создаём длинное вещественное число
	cout << "a = " << a << endl; // выводим созданное число
}

Добавляем арифметику

Вряд ли создавать и выводить числа будет долгое время интересным, поэтому переходим к, пожалуй, самой важной части - добавлению арифметических операций, а именно:

• сложение (LongDouble operator+(const LongDouble& x) const)
• вычитание (LongDouble operator-(const LongDouble& x) const)
• умножение (LongDouble operator*(const LongDouble& x) const)
• деление (LongDouble operator/(const LongDouble& x) const)

Умножение двух длинных вещественных чисел

Начнём с умножения, так как это наиболее простая операция для вещественных чисел. Чтобы умножить два вещественных числа, нужно умножить два целых числа, записанных в векторе с помощью умножения в столбик, а затем сложить экспоненты и перемножить знаки. Да, вот так просто!

LongDouble LongDouble::operator*(const LongDouble& x) const {
	size_t len = digits.size() + x.digits.size(); // максимальная длина нового числа не больше суммы длин перемножаемых чисел

	LongDouble res; // создадим объект для результата

	res.sign = sign * x.sign; // перемножаем знаки
	res.digits = vector<int>(len, 0); // создаём вектор из нулей
	res.exponent = exponent + x.exponent; // складываем экспоненты

	// умножаем числа в столбик
	for (size_t i = 0; i < digits.size(); i++)
		for (size_t j = 0; j < x.digits.size(); j++)
			res.digits[i + j + 1] += digits[i] * x.digits[j]; // прибавляем результат произведения цифр из i-го и j-го разрядов

	// в результате такого перемножения в ячейках могли получиться двузначные числа, поэтому нужно выполнить переносы
	for (size_t i = len - 1; i > 0; i--) {
		res.digits[i - 1] += res.digits[i] / 10; // добавляем к более старшему разряду десятки текущего разряда
		res.digits[i] %= 10; // оставляем только единицы у текущего разряда
	}

	res.removeZeroes(); // удаляем нули, как в конструкторе из строки, если таковые имеются

	return res; // возвращаем результат умножения двух чисел
}

Сложение и вычитание длинных вещественных чисел

Как ни странно, но сложение и вычитание двух вещественных чисел будет сложнее умножения, так как нужно дополнительно обрабатывать числа с разными знаками, отправляя на противоположную операцию, а также выравнивать числа по разрядам.

Но перед тем, как реализовывать сложение и вычитание, добавим в наш класс проверку на то, что одно число больше другого и получение числа с обратным знаком (унарный минус):

// унарный минус
LongDouble LongDouble::operator-() const {
	LongDouble res; // создаём число
	res.sign = -sign; // меняем знак на противоположный
	res.exponent = exponent; // копируем экспоненту
	res.digits = vector<int>(digits); // копируем цифры

	return res; // возвращаем результат
}

// проверка, что число больше другого
bool LongDouble::operator>(const LongDouble& x) const {
	if (sign != x.sign)
		return sign > x.sign; // если знаки разные, то положительное число больше

	if (exponent != x.exponent)
		return (exponent > x.exponent) ^ (sign == -1); // если экспоненты разные, то больше число с большей экспонентой с учётом знака
	
	// копируем вектора
	vector<int> d1(digits);
	vector<int> d2(x.digits);
	size_t size = max(d1.size(), d2.size()); // определяем максимальный размер векторов
	
	// выравниваем размеры векторов, добавляя в концы нули
	while (d1.size() != size)
		d1.push_back(0);

	while (d2.size() != size)
		d2.push_back(0);
	
	// проходим по всем цифрам числа
	for (size_t i = 0; i < size; i++)
		if (d1[i] != d2[i])
			return (d1[i] > d2[i]) ^ (sign == -1); // если в каком-то разряде цифры разные, то больше число с большей цифрой с учётом знака

	return false; // иначе числа равны, а значит не больше
}

А вот теперь действительно сложение и вычитание

Теперь мы действительно готовы реализовать сложение и вычитание. Начнём со сложения. Алгоритм будет таков: если знаки совпадают, то выравниваем разряды чисел путём задания обоим числам экспоненты, равной максимальной из них, дополняя другое число нулями в начале, складываем полученные цифры, проверяем возможные переполнения и избавляемся от нулей. Если же знаки разные, то если первое число отрицательно, то выполняем вычитание из второго числа первого, умноженного на -1, а иначе вычитаем из первого второе, умноженное на -1.

LongDouble LongDouble::operator+(const LongDouble& x) const {
	if (sign == x.sign) { // если знаки одинаковые
		long exp1 = exponent; // запоминаем экспоненту первого
		long exp2 = x.exponent; // и второго чисел
		long exp = max(exp1, exp2); // находим максимальную экспоненту. К ней будем приводить оба числа

		// создаём вектора из векторов цифр наших чисел
		vector<int> d1(digits);
		vector<int> d2(x.digits);
		
		// выравниваем экспоненты
		while (exp1 != exp) {
			d1.insert(d1.begin(), 0); // добавляя нули в начало первого
			exp1++;
		}

		while (exp2 != exp) {
			d2.insert(d2.begin(), 0); // и в начало второго векторов
			exp2++;
		}

		size_t size = max(d1.size(), d2.size()); // определяем максимальный размер векторов

		// выравниваем размеры векторов, добавляя нули в конец каждого из них
		while (d1.size() != size)
			d1.push_back(0);

		while (d2.size() != size)
			d2.push_back(0);

		size_t len = 1 + size;

		LongDouble res; // создаём новое число

		res.sign = sign; // знак результата совпадает со знаком чисел
		res.digits = vector<int>(len, 0); // создаём вектор из len нулей

		// заполняем каждую ячейку вектора суммой двух соответствующих цифр
		for (size_t i = 0; i < size; i++)
			res.digits[i + 1] = d1[i] + d2[i];

		// проверяем переполнения
		for (size_t i = len - 1; i > 0; i--) {
			res.digits[i - 1] += res.digits[i] / 10; // прибавляем к более старшему разряду десятки текущего
			res.digits[i] %= 10; // оставляем у текущего разряда только единицы
		}

		res.exponent = exp + 1; // восстанавливаем экспоненту
		res.removeZeroes(); // убираем нули

		return res; // возвращаем число
	}

	if (sign == -1)
		return x - (-(*this)); // если первое число отрицательное, то из второго вычитаем первое с обратным знаком
	
	return *this - (-x); // иначе из первого вычитаем второе с обратным знаком
}

В случае с вычитанием, нам потребуется узнать, какое из чисел больше, чтобы из большего числа вычитать меньшее, однако в остальном логика будет такая же, как и в случае со сложением.

LongDouble LongDouble::operator-(const LongDouble& x) const {
	if (sign == 1 && x.sign == 1) { // если боа числа положительны
		bool cmp = *this > x; // получаем флаг того, больше ли первое число

		long exp1 = cmp ? exponent : x.exponent; // сохраняем экспоненту большего числа
		long exp2 = cmp ? x.exponent : exponent; // сохраняем экспоненту меньшего числа
		long exp = max(exp1, exp2); // определяем максимальную экспоненту, чтобы к ней привести числа

		vector<int> d1(cmp ? digits : x.digits); // запоминаем вектор цифр большего числа
		vector<int> d2(cmp ? x.digits : digits); // запоминаем вектор цифр меньшего числа

		// выравниваем экспоненты как при сложении (добавляя нули вначале числа)
		while (exp1 != exp) {
			d1.insert(d1.begin(), 0);
			exp1++;
		}

		while (exp2 != exp) {
			d2.insert(d2.begin(), 0);
			exp2++;
		}

		size_t size = max(d1.size(), d2.size()); // определяем максимальный размер

		// добавляем нули в конец векторов цифр
		while (d1.size() != size)
			d1.push_back(0);

		while (d2.size() != size)
			d2.push_back(0);

		size_t len = 1 + size;

		LongDouble res; // создаём число для результата

		res.sign = cmp ? 1 : -1; // знак будет 1, если первое больше второго, и -1, если первое меньше второго
		res.digits = vector<int>(len, 0); // создаём вектор из len нулей

		for (size_t i = 0; i < size; i++)
			res.digits[i + 1] = d1[i] - d2[i]; // вычитаем соответствующие разряды

		// обрабатываем переполнения
		for (size_t i = len - 1; i > 0; i--) {
			if (res.digits[i] < 0) { // если текущий разряд стал меньше нуля
				res.digits[i] += 10; // занимаем у предыдущего, прибавляя 10 к текущему
				res.digits[i - 1]--; // уменьшаем на 1 предыдущий разряд
			}
		}

		res.exponent = exp + 1; // восстанавливаем экспоненту
		res.removeZeroes(); // удаляем лишнии нули

		return res; // возвращаем результат
	}

	if (sign == -1 && x.sign == -1)
		return (-x) - (-(*this)); // если оба числа отрицательны, то из второго с обратным знаком вычитаем первое с обратным знаком
	
	return *this + (-x); // если знаки разные, прибавляем к первому отрицательное второе
}

Поздравляю! Теперь мы можем складывать, вычитать и умножать вещественные числа! Это гораздо больше, чем только создавать и выводить, не правда ли? Но для полного счастья ещё хочется реализовать деление чисел. Что ж, давайте реализовывать!

Деление двух длинных вещественных чисел

Что бы разделить одно вещественное число на другое, будем использовать следующую стратегию: найдём обратный элемент (1 / x) к числу, на которое хотим поделить и умножим его на число, которое хотели разделить. То есть нам дополнительно потребуется функция поиска обратного значения, которую логично сделать даже публичной, чтобы пользователь всегда мог найти обратный элемент самостоятельно, не прибегая к делению 1 на число.

При нахождении обратного числа, может получиться ситуация, что разделить до конца не получится (например, 1 / 3), поэтому ограничим максимальное количество цифр после запятой для нахождения обратного элемента, заведя константу divDigits и положим её значение равным, например, 1000. Весьма неплохая такая точность! :)

LongDouble LongDouble::inverse() const {
	if (digits.size() == 1 && digits[0] == 0)
		throw string("LongDouble LongDouble::inverse() - division by zero!"); // делить на ноль нельзя, поэтому бросим исключение

	LongDouble x(*this); // скопируем число,
	x.sign = 1; // сделав его положительным

	LongDouble d("1"); // создадим то, что будем делить на x

	LongDouble res; // создадит объект для результата
	res.sign = sign; // знак результата совпадёт со знаком числа
	res.exponent = 1; // начнём с единичной экспоненты
	res.digits = vector<int>(); // создадим пустой вектор для цифр обратного элемента

	// пока число меньше 1
	while (x < 1) {
		x.exponent++; // будем увеличивать его экспоненту (умножать на 10 фактически)
		res.exponent++; // и заодно экспоненту результата
	}
	
	// дальше сдлеаем число d большим x, также умножая его на 10, чтобы получить число 100...0
	while (d < x)
		d.exponent++;	
	
	res.exponent -= d.exponent - 1; // подсчитаем реальное количество цифр в целой части

	size_t numbers = 0; // количество уже вычисленных цифр дробной части
	size_t totalNumbers = divDigits + max((long) 0, res.exponent); // количество цифр с учётом целой части

	do {
		int div = 0; // будущая цифра
		
		// считаем, сколько раз нужно вычесть x из d
		while (d >= x) {
			div++; 
			d = d - x;
		}

		// увеличиваем остаток в 10 раз
		d.exponent++; 
		d.removeZeroes();
		res.digits.push_back(div); // записываем сформированную цифру
		numbers++; // увеличиваем число вычисленных цифр
	} while (d != 0 && numbers < totalNumbers); // считаем до тех пор, пока не дойдём до нулевого остатка или пока не превысим точность

	return res; // возвращаем результат
}

А теперь можно и деление писать:

LongDouble LongDouble::operator/(const LongDouble& x) const {
	LongDouble res = *this * x.inverse(); // выполняем умножение на обратный элемент

	// избавляемся от записи вида d...dd(9), которая эквивалентна d...dd+1
	size_t i = res.digits.size() - 1 - max((long)0, exponent);
	size_t n = max((long) 0, res.exponent);

	// если в указанном месте девятка, то ищем место, в котором девятки закончатся
	if (i > n && res.digits[i] == 9) {
		while (i > n && res.digits[i] == 9)
			i--;
		
		// если дошли до целой части
		if (res.digits[i] == 9) {
			res.digits.erase(res.digits.begin() + n, res.digits.end()); // то удаляем всю вещественную часть
			res = res + res.sign; // и прибавляем 1 (или -1 к отрицательному) 
		}
		else {
			res.digits.erase(res.digits.begin() + i + 1, res.digits.end()); // иначе обрезаем с найденного места
			res.digits[i]++; // и увеличиваем на 1 текущий разряд
		}
	}

	return res; // возвращаем результат деления
}

Теперь у нас есть настоящая длинная вещественная арифметика, которую можно складывать, вычитать, умножать и делить! Ну и выводить, конечно же! Чего теперь не хватает для полного счастья? Пожалуй, методов для сравнения чисел. Метод "больше" уже создан, а значит осталось реализовать только 5 методов сравнения: "меньше", "больше или равно", "меньше или равно" и "не равно". Значит реализуем!

Методы для сравнения длинных вещественных чисел

bool LongDouble::operator==(const LongDouble& x) const {
	if (sign != x.sign)
		return false;

	if (exponent != x.exponent)
		return false;

	if (digits.size() != x.digits.size())
		return false;

	for (size_t i = 0; i < digits.size(); i++)
		if (digits[i] != x.digits[i])
			return false;

	return true;
}

bool LongDouble::operator<(const LongDouble& x) const {
	return !(*this > x || *this == x);
}

bool LongDouble::operator>=(const LongDouble& x) const {
	return (*this > x || *this == x);
}

bool LongDouble::operator<=(const LongDouble& x) const {
	return *this < x || *this == x;
}

bool LongDouble::operator!=(const LongDouble& x) const {
	return !(*this == x);
}

Пожалуй, на этом остановимся, ведь по сути класс завершён! Для тех, кому хочется самостоятельно протестировать результаты наших трудов, переходите по ссылке на гитхаб страницу проекта. Там вы сможете найти файлы с классом, а также некоторые методы, которые не были описаны в этой статье. К тому же, в проекте находится основная программа, выполняющая тестирование всех арифметических методов с замером времени выполнения операций!