/*
 *  Copyright 2007 Daniel Blazewicz
 *
 *  e-mail: klajok @ interia . pl
 */

/*
 *      This program is free software; you can redistribute it and/or modify
 *      it under the terms of the GNU General Public License as published by
 *      the Free Software Foundation; either version 2 of the License, or
 *      (at your option) any later version.
 *
 *      This program is distributed in the hope that it will be useful,
 *      but WITHOUT ANY WARRANTY; without even the implied warranty of
 *      MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *      GNU General Public License for more details.
 *
 *      You should have received a copy of the GNU General Public License
 *      along with this program; if not, write to the Free Software
 *      Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
 *
 *      Niniejszy program jest wolnym oprogramowaniem; możesz go
 *      rozprowadzać dalej i/lub modyfikować na warunkach Powszechnej
 *      Licencji Publicznej GNU, wydanej przez Fundację Wolnego
 *      Oprogramowania - według wersji 2 tej Licencji lub (według twojego
 *      wyboru) którejś z późniejszych wersji.
 *
 *      Niniejszy program rozpowszechniany jest z nadzieją, iż będzie on
 *      użyteczny - jednak BEZ JAKIEJKOLWIEK GWARANCJI, nawet domyślnej
 *      gwarancji PRZYDATNOŚCI HANDLOWEJ albo PRZYDATNOŚCI DO OKREŚLONYCH
 *      ZASTOSOWAŃ. W celu uzyskania bliższych informacji sięgnij do
 *      Powszechnej Licencji Publicznej GNU.
 *
 *      Z pewnością wraz z niniejszym programem otrzymałeś też egzemplarz
 *      Powszechnej Licencji Publicznej GNU (GNU General Public License);
 *      jeśli nie - napisz do Free Software Foundation, Inc., 59 Temple
 *      Place, Fifth Floor, Boston, MA  02110-1301  USA
 */

#include <stdlib.h>
#include <iostream>
#include <gmp.h>
#include <gmpxx.h>
#include <time.h>

int main(int argc, char* argv[]){

  /*--------------------Dane wejsciowe---------------------*/

  if(argc < 2){
    std::cerr << "Skladnia wywolania:\n";
    std::cerr << "   izomery maksymalna_liczba_atomow_C\n";
    std::cerr << "   np.  izomery 100\n";
    std::cerr << "lub\n";
    std::cerr << "   izomery zakres_liczby_atomow_C\n";
    std::cerr << "   np.  izomery 90 100\n";
    return 1;
  }
  long N = atol(argv[1]);
  long M = 1;
  if (argc > 2)
    M = atol(argv[2]);
  if (M > N) {
    long temp = M;
    M = N;
    N = temp;
  }
  if (M < 1) {
  	std::cerr << "Podałeś co najmniej jeden nieprawidłowy parametr wywołania!\n";
  	return 1;
  }
  if (N > 400000) {
    std::cerr << "Chyba oszalalas/oszalales!!\n";
    return 0;
  }

  
  /*------------Obliczenia dla grup alkilowych------------*/

  mpz_class* alkil = new mpz_class[N/2+1];
  mpz_class* alkil2 = new mpz_class[N/2+1];
  mpz_class* jpom = new mpz_class[N/2+1];

  alkil[0] = 1;
  alkil2[0] = 1;
  jpom[0] = 0;
  for(long i=1; i<=N/2; ++i){
    alkil[i] = 0;
    jpom[i] = 0;
  }

  std::cerr << "Zaczynam liczyc...\n";
  clock_t czas = std::clock();

  std::cout << "Liczba izomerow grup alkilowych CnH2n+1 :\n";

  for(long n=1; n<=N/2; ++n){
    
    for (long k = (n+4)/3; k < n-1; ++k)
      alkil[n] += jpom[n-k-1] * alkil[k];
    
    for (long k = 0, a = n-1; a >= 0; ++k, a-=2)
      if (k != a)
	alkil[n] += alkil2[k] * alkil[a];
      else
	alkil[n] += (alkil2[k] * (alkil[k]+2)) / 3;
    
    alkil2[n] = ((alkil[n]+1) * alkil[n]) / 2;

    for (long k = n/2; k < (2*n)/3; ++k)
      jpom[k] += alkil[2*k-n+1] * alkil[n-k-1];
    
    std::cout << n << "  " << alkil[n] << "\n";
  }

  std::cerr << "Czas obliczen dla grup alkilowych: " << (std::clock()-czas) / 1e6 << " s.\n";
  czas = std::clock();

  /*---------------Obliczenia dla alkanow---------------*/

  std::cout << "Liczba izomerow alkanow CnH2n+2 :\n";
  
  mpz_class * alkan = new mpz_class[N+1];

  for(long n=M; n<=N; ++n)
    if (n % 2 == 0)
      alkan[n] = alkil2[n/2];
    else
      alkan[n] = 0;

  for (long k = 1; k <= 2*((N-1)/2)-1; ++k){

    mpz_class jk = 0;

    long kmin = std::max((M-k+3)/2, std::max((k+6)/4, k-(N-1)/2));
    long k1r = ((k%2) == 0) ? k-M/2 : k-(M-1)/2;
    long k1l = ((k%2) == 0) ? k-(N-1)/2 : k-(N-2)/2;
    for (long k1 = (k-1)/2; k1 >= kmin; --k1){

      jk += alkil[k1] * alkil[k-k1];

      long kpom = 2*(k-k1)+1;

      for (long n = std::max(kpom, M); n <= std::min(N, kpom+1); ++n)
	if (k1 < n-k){
	  mpz_class suma_b = 0;
	  for (long b = k1-1; b >= (n-k+1)/2; --b)
	    suma_b += alkil[n-k-1-b] * alkil[b];
	  alkan[n] += suma_b * jk;
	} else
	  alkan[n] += jpom[n-k-1] * jk;
	  
      for (long n = std::max(k+k1,std::max(kpom+2, M)); n <= std::min(N, k+2*k1-2); ++n)
	alkan[n] += alkil[n-k-k1] * alkil[k1-1] * jk;

      if ((k1l <= k1) && (k1 <= k1r)) {
	long n = k+2*((k+1)/2)+1-2*k1;
	alkan[n] += alkil2[(n-k-1)/2] * jk;
      }

    }

    kmin = std::min(kmin, (k+1)/2);
    if (k1l <= kmin-1) {
      if (k1r < kmin-1) {
	for (long k1 = kmin-1; k1 > std::max(k1r, long(-1)); --k1)
	  jk += alkil[k1] * alkil[k-k1];
	kmin = k1r+1;
      }
      for (long k1 = kmin-1; k1 >= k1l; --k1) {
	if (k1 >= 0)
	  jk += alkil[k1] * alkil[k-k1];
	long n = k+2*((k+1)/2)+1-2*k1;
	alkan[n] += alkil2[(n-k-1)/2] * jk;
      }
    }

  }
  
  std::cerr << "Czas glownych obliczen dla alkanow: " << (std::clock()-czas) / 1e6 << " s.\n";
  czas = std::clock();

  for(long n=M; n<=N; ++n) {

    for (long a = (n-1)/3, b = n-1-3*a; b <= (n-1)/2; --a, b+=3)
      if (b != a)
	alkan[n] -= ((alkil2[a] * (2*alkil[a]-2))/3) * alkil[b];

    if ( (n-1) % 2 == 0 )
      for (long a = 0,  b = (n-1)/2-a; b > a; ++a, --b)
	alkan[n] += alkil2[a] * alkil2[b];

    if ( (n-1) % 4 == 0 ) {
      long a = (n-1)/4;
      alkan[n] += ((alkil2[a] + alkil[a]) * (alkil2[a] + alkil[a] + 1)) / 6;
    }

    std::cout << n << "  " << alkan[n] << "\n";

  }

  std::cerr << "Czas drobniejszych obliczen dla alkanow: " << (std::clock()-czas) / 1e6 << " s.\n";

  delete [] alkan;
  delete [] jpom;
  delete [] alkil;
  delete [] alkil2;

  return 0;
}
