// Algorithms from Theorem 2
// including the prefix sums AND the memory optimization
// Using the GMPlib big integers from the CGAL library
// Runtime O(n^4), Memory O(n^3) times the bigint-overhead (roughly another O(nlogn)).

// Just compile and run it to get a quick demo or look at the main function for some examples:
// Compilation steps with CGAL:
// Run: cgal_create_cmake_script
// To the file CMakeLists.txt add this line: set(CMAKE_CXX_FLAGS ${CMAKE_CXX_FLAGS} " -std=c++11")
// Run: cmake .
// Run: make
// Execute: ./motzkin_leftright_bignum


#include <iostream>
#include <vector>
#include <map>
#include <algorithm>
#include <cassert>
#include <CGAL/basic.h>
#include <CGAL/Gmpz.h>
typedef CGAL::Gmpz bignum;

using namespace std;
using in = long long int;
using PII = pair<in,in>;
using VI = vector<in>;
using VB = vector<bignum>;

#define sz(v) (v).size()
#define forn(i,n) for(in i=0; i<(n); ++i)
#define forv(i,v) forn(i,sz(v))
#define DEB(x)  // Replace by "DEB(x) x" to get debug output
using namespace std;

bignum tob(in x) {
	return (bignum)(double)x;
}

// M(n,A): Unweighted Version
vector<map<VI, bignum>> LFMemo;
vector<map<VI, bignum>> PSLFMemo; // Prefix sums
// D(n,d): Weighted Version
vector<map<VI, bignum>> WLFMemo;
vector<map<VI, bignum>> PSWLFMemo; // Prefix sums

int globalw;

const in MAXW = 110;
void init() {
	LFMemo.resize(MAXW);
	PSLFMemo.resize(MAXW);
	WLFMemo.resize(MAXW);
	PSWLFMemo.resize(MAXW);
}

bignum LFcount(in A, in w, in l);

bignum SumLFcount(in A, in w, in la, in lb) { // sum for l in [la,lb]
	if(la > lb) return 0;
	if(la < 0) return 0;
	if(la==lb) return LFcount(A,w,la);
	if(PSLFMemo[w].find({A,w,lb}) == PSLFMemo[w].end()) {
		PSLFMemo[w][{A,w,lb}] = SumLFcount(A,w,0,lb-1) + LFcount(A,w,lb);
	}
	return PSLFMemo[w][{A,w,lb}] - SumLFcount(A,w,0,la-1);
}

bignum LFcount(in A, in w, in l) {
	if(w>3) {
		LFMemo[w-3].clear();
		PSLFMemo[w-3].clear();
	}
	if(A==0 && w == 0 && l==0 ) {
		return 1;
	}
	if(A < 0 || w <= 0 || l < 0) return 0;
	if(A > w*w/2) return 0;
	if(LFMemo[w].find({A,w,l}) != LFMemo[w].end()) {
		return LFMemo[w][{A,w,l}];
	}
	if(w + 5 < globalw)
		cout << "A=" << A << " w=" << w << " l=" << l << endl;
	// if(w>5) {
	// 	// cout << "del for w " << w << endl;
	// 	PSLFMemo[w-5].clear();
	// 	LFMemo[w-5].clear();
	// }
	bignum res = 0;
	res += SumLFcount(A-l, w-1, l, w/2);
	// for(in ll = l; ll<=w/2; ll++) {
	// 	res += LFcount(A-l, w-1, ll);
	// }
	if(l>0) {
		res += SumLFcount(A-(2*l-1), w-2, l-1, w/2);
	    //	for(in ll = l-1; ll<=w/2; ll++) {
		// 	res += LFcount(A-(2*l-1), w-2, ll);
		// }
	}
	LFMemo[w][{A,w,l}] = res;
	// cout << " computed " << "LFcount(" << A << "," << w << "," << l << ") = " << res << endl;
	return res;
}
bignum LFcount(in A, in w) {
	bignum res = 0;
	for(in l = 0; l <= w; l++) {
		bignum x =  LFcount(A,w,l);
		// cout << "add " << "LFcount(" << A << "," << w << "," << l << ") =... " << x <<  endl;
		res += x;
	}
	return res;
}
bignum LFcount(in w) {
	bignum res = 0;
	for(in A = 0; A <= w*w; A++) {
		res += LFcount(A,w);
	}
	return res;
}

// Weighted Version
bignum WLFcount(in A, in w, in l);

bignum SumWLFcount(in A, in w, in la, in lb) { // sum for l in [la,lb]
	if(la > lb) return 0;
	if(la < 0) return 0;
	if(la==lb) return WLFcount(A,w,la);
	if(PSWLFMemo[w].find({A,w,lb}) == PSWLFMemo[w].end()) {
		PSWLFMemo[w][{A,w,lb}] = SumWLFcount(A,w,0,lb-1) + WLFcount(A,w,lb);
	}
	return PSWLFMemo[w][{A,w,lb}] - SumWLFcount(A,w,0,la-1);
}

bignum WLFcount(in A, in w, in l) {
	if(w>3) {
		LFMemo[w-3].clear();
		PSLFMemo[w-3].clear();
	}
	// cout << "A " << A << " w " << w << " l " << l << endl;
	if(A==0 && w == 0 && l==0 ) {
		// cout << "base" << endl;
		return 1;
	}
	if(A < 0 || w <= 0 || l < 0) return 0;
	if(A > w*w/2) return 0;
	if(WLFMemo[w].find({A,w,l}) != WLFMemo[w].end()) {
		return WLFMemo[w][{A,w,l}];
	}
	bignum res = 0;
	res += tob(2*l+1)*SumWLFcount(A-l, w-1, l, w/2);
	if(l>0) {
		res += tob(l)*tob(l)*SumWLFcount(A-(2*l-1), w-2, l-1, w/2);
	}
	WLFMemo[w][{A,w,l}] = res;
	return res;
}
bignum WLFcount(in A, in w) {
	bignum res = 0;
	for(in l = 0; l <= w; l++) {
		bignum x =  WLFcount(A,w,l);
		res += x;
	}
	return res;
}
bignum WLFcount(in w) {
	bignum res = 0;
	for(in A = 0; A <= w*w; A++) {
		res += WLFcount(A,w);
	}
	return res;
}

int main() {
	init();

	in A,w,h;
	in WStop = 100;

	cout << "\nMotzkin Path Counting Demo\n\n";

	cout << "Motzkin Numbers from 1 to 10\n";
	cout << "see: https://en.wikipedia.org/wiki/Motzkin_number\n";
	cout << "see: http://oeis.org/A001006\n";

	for(in w = 1; w<=WStop; w++) {
		bignum res = LFcount(w);
		cout << "M(" << w << ") = " << res << endl;
	}

	cout << "\n\nM(n,A): Print Motzkin Path Counts by width and area using the last fall DP\n";
	cout << "see: http://oeis.org/A129181\n";
	for(w = 1; w<=WStop; w++) {
		for(A = 0; A <= w*w; A++) {
			bignum res = LFcount(A,w);
			// assert(res == weighted_count(A,w));
			if(res>0) cout << res << "\t";
		}
		cout << endl;
	}

	cout << "\n\nPrint Factorial Numbers w! = Sum of Weighted Motzkin Path Counts by width for a fixed w\n";
	cout << "see: http://oeis.org/A000142\n";
	for(w = 1; w<=100; w++) {
		globalw = w;
		bignum res = WLFcount(w);
		cout << w << "! = " << res << endl;
	}

	cout << "\n\nD(n,A): Print Weighted Motzkin Path Counts by width and area using the last fall DP\n";
	cout << "see: https://oeis.org/A062869\n";
	for(w = 1; w<=WStop; w++) {
		for(A = 0; A <= w*w; A++) {
			bignum res = WLFcount(A,w);
			// assert(res == weighted_count(A,w));
			if(res>0) cout << res << "\t";
		}
		cout << endl;
	}
}