/*************************************************************************** * Copyright (c) 2017 * Joel Berman & Peter Koehler * * Project: FDn * Filename: a132581.cpp * Purpose: Caclulation Algorithm for A132581 * Author: Peter Koehler * * - HoC - * Date | Author | Description * --------+--------+-------------------------------- * 05.07.17 pk Submission to OEIS ***************************************************************************/ //! \file fdseqownmpi.cpp //! //! This program uses a straightforward implementation of the "Splitting result" of //! J.Berman and P.Koehler, Cardinalities of finite distributive lattices, //! Mitteilungen aus dem Mathematischen Seminar Giessen, 121 (1976), 103 - 124. //! which is implemented in the function f below. //! //! The function f is called for the numbers 2^i-1, i=0,...,n, where n is the last //! index of the sequence to be computed. The numbers 2^i-1 correspond to the posets //! of the sequence. If i is of the form 2^k for some k, then the poset corresponding //! to 2^i-1 is the Boolean algebra with i elements. The computed values are printed //! to Standard Output. //! //! As f is called recursively, for speedup the computed values are stored in a map, //! where they can be looked up efficiently if used in subsequent steps. //! The number of entries of the map is also printed to Standard Output after the last //! step. //! //! The Multi-precision integer arithmetic is implemented in a separate class //! BigInt - it would also be possible to use the GNU multiprecion package GMP. #define _CRT_SECURE_NO_WARNINGS #include <iostream> #include <map> #include <vector> #include <cstdio> #include <cstdlib> #include <cstring> using namespace std; //! \class BigInt //! //! The class BigInt provides a multiprecision //! integer arithmetic. Oppeosed to e.g. GNU //! GMP it uses a (compile time) fixed bitsize, in order //! to make the operations simple and fast. class BigInt { static const int Dim = 16; // 256 Bit - BigEndian public: BigInt(unsigned int n = 0) { memset(nums, 0, sizeof(nums)); nums[Dim - 1] = n & 0xffff; nums[Dim - 2] = (n >> 16) & 0xffff; } BigInt(const BigInt & other) { memcpy(nums, other.nums, sizeof(nums)); } BigInt operator+ (const BigInt & other); BigInt & operator+= (const BigInt & other); bool operator!= (int val); BigInt operator<< (unsigned int count); BigInt & operator>>= (unsigned int count); BigInt operator& (const BigInt & other); BigInt & operator&= (const BigInt & other); BigInt operator ~(); BigInt operator/ (unsigned short value); unsigned short operator% (unsigned short value); BigInt & andComp(const BigInt & other); private: unsigned short nums[Dim]; static void div(const BigInt & first, unsigned short second, BigInt "ient, unsigned short & remainder); static void print(std::ostream & stream, BigInt arg); friend bool operator< (const BigInt & first, const BigInt & second); friend std::ostream & operator<<(std::ostream &, BigInt); }; BigInt BigInt::operator+ (const BigInt & other) { BigInt result; unsigned short carry = 0; for (int n = Dim - 1; n >= 0; n--) { unsigned int sum = carry + nums[n] + other.nums[n]; result.nums[n] = sum & 0xffff; carry = sum >= 0x10000 ? 1 : 0; } return result; } BigInt & BigInt::operator+= (const BigInt & other) { unsigned short carry = 0; for (int n = Dim - 1; n >= 0; n--) { unsigned int sum = carry + nums[n] + other.nums[n]; nums[n] = sum & 0xffff; carry = sum >= 0x10000 ? 1 : 0; } return *this; } bool BigInt::operator!= (int val) { for (int n = 0; n < Dim - 2; n++) { if (nums[n] != 0) return true; } return (nums[Dim - 1] != (val & 0xffff)) || (nums[Dim - 2] != ((val >> 16) & 0xffff)); } BigInt BigInt::operator<< (unsigned int count) { BigInt result; int offset = count / 16; unsigned int shift = count % 16; for (int n = 0; n < Dim; n++) { unsigned shiftValue = nums[n] << shift; if (n - offset >= 0) { result.nums[n - offset] |= shiftValue & 0xffff; if (n - offset > 0) result.nums[n - offset - 1] |= (shiftValue >> 16) & 0xffff; } } return result; } BigInt & BigInt::operator>>= (unsigned int count) { BigInt temp; unsigned int offset = count / 16; unsigned int shift = count % 16; for (int n = 0; n < Dim; n++) { unsigned shiftValue = (((unsigned)nums[n]) << 16) >> shift; if (n + offset < Dim) { temp.nums[n + offset] |= (shiftValue >> 16) & 0xffff; if (n + offset < Dim - 1) temp.nums[n + offset + 1] |= (shiftValue & 0xffff); } } *this = temp; return *this; } BigInt BigInt::operator& (const BigInt & other) { BigInt result; for (int n = 0; n < Dim; n++) result.nums[n] = nums[n] & other.nums[n]; return result; } BigInt & BigInt::operator&= (const BigInt & other) { for (int n = 0; n < Dim; n++) nums[n] &= other.nums[n]; return *this; } BigInt BigInt::operator ~() { BigInt result; for (int n = 0; n < Dim; n++) result.nums[n] = (~nums[n]) & 0xffff; return result; } bool operator< (const BigInt & first, const BigInt & second) { int index = 0; while (first.nums[index] == second.nums[index] && index < BigInt::Dim) index++; return index != BigInt::Dim && first.nums[index] < second.nums[index]; } BigInt BigInt::operator/ (unsigned short value) { unsigned short remainder; BigInt result; div(*this, value, result, remainder); return result; } unsigned short BigInt::operator% (unsigned short value) { unsigned short remainder; BigInt result; div(*this, value, result, remainder); return remainder; } //! \brief And-With-Complement (x &= ~y) BigInt & BigInt::andComp(const BigInt & other) { for (int n = 0; n < Dim; n++) nums[n] &= ~other.nums[n]; return *this; } void BigInt::div(const BigInt & first, const unsigned short second, BigInt & quotient, unsigned short & remainder) { unsigned int carry = 0; for (int n = 0; n < Dim; n++) { unsigned int temp = 0x10000 * carry + first.nums[n]; quotient.nums[n] = temp / second; carry = temp % second; } remainder = carry & 0xffff; } void BigInt::print(ostream & stream, BigInt arg) { if (arg < 10) stream << arg.nums[Dim - 1]; else { unsigned short remainder; div(arg, 10, arg, remainder); stream << arg << remainder; } } std::ostream & operator<<(std::ostream & stream, BigInt arg) { BigInt::print(stream, arg); return stream; } typedef std::map<BigInt, BigInt> ValueMap; static ValueMap valueMap; static vector<BigInt> powers; static unsigned int bitsize; //! \brief Computation Function //! \param m: Poset expressed as number //! \result: number of semiideals of the poset BigInt f ( BigInt m ) { register unsigned x, t ; BigInt u, s, p; // If the value is known, take it from the map ValueMap::iterator iter = valueMap.find(m); if ( iter != valueMap.end() ) { return iter->second; } u = m; // Determine log2(m) for (x = 0; u != 0; x++) { u >>= 1; } x -= 1; // Size of Poset P-{x} BigInt p1 = m; p1.andComp(powers[x]); s = f(p1); // Size of Poset P-cone(x) p = m; for ( t = 0; t <= x; t++ ) { if ( (t & x) == t ) { p.andComp(powers[t]); } } s += f(p); // Store the value valueMap[m] = s; return s; } //! \brief Main Function //! \param argc: Number of command line arguments //! \param args: Argument strings int main (int argc, char *argv[] ) { int x = -1 ; if ( ( argc == 1 ) || ( *argv[1] == '?' ) || ( 1 != sscanf ( argv [1], "%d", &x ) ) || ( x < 0 ) ) cerr << "The program calculates the first terms of\n" "the integer sequence A132581 of OEIS" << endl; else { bitsize = x; BigInt n(1); int i; valueMap[0] = BigInt(1); valueMap[1] = BigInt(2); for(int n=0; n<x; n++) { powers.push_back(BigInt(1) << n); } cout << "# First " << x << " terms of A132581 - computed by Peter Koehler (05.07.2017)" << endl; for (i = 1; i <= x; i = i + 1, n = n + n + BigInt(1)) { BigInt cur = f(n); std::cout << i << " " << cur << std::endl; } cout << "# Map usage: " << valueMap.size() << " entries" << endl;; } return 0; }