/*
 This C++ program counts and displays the
 1 1 2 5 15 54 212 908 4011 18260
 incongruental non-labbeled undirected trees on square
 lattices with n=1,2,3,.., nodes and n-1 edges.
 http://oeis.org/A056841

 Richard J. Mathar, 2006-04-09
*/
#include <limits.h>
#include <iostream>
#include <vector>
#include <algorithm>
#include <unistd.h>

using namespace std ;

/** A node is one point on the square lattice.
* A node contains simply the two integer coordinates of a point on the
* square lattice.
*/
class node {
public:
/	/* xy[0] is the x, and xy[1] the y-component of the coordinates */
	short xy[2] ;

	node()
	{
		xy[0]=xy[1] = SHRT_MIN ;
	}

	/** Constructor with a pair of coordinates.
	* @param x the x coordinate of the point
	* @param y the y coordinate of the point
	*/
	node(short x, short y)
	{
		xy[0]=x ;
		xy[1] = y ;
	}

	/** Constructor with a given original node translated into a direction.
	* @param n the node to start from
	* @param dir a parameter from 0 to 4 to indicate in which
	*   direction of the node \c n this node that is created here
	*   is to be placed.
	*/
	node(const node & n, int dir)
	{
		/** For values of the \c dir argument from 0 to 3
		* we create the new node one step into the East (+x),
		* North (+y), West (-x) or South (-y) from the original one.
		*/
		switch(dir)
		{
		case 0:
			xy[0] = n.xy[0]+1;
			xy[1] = n.xy[1] ;
			break ;
		case 1:
			xy[0] = n.xy[0];
			xy[1] = n.xy[1]+1 ;
			break ;
		case 2:
			xy[0] = n.xy[0]-1;
			xy[1] = n.xy[1] ;
			break ;
		case 3:
			xy[0] = n.xy[0];
			xy[1] = n.xy[1]-1 ;
			break ;
		}
	}

	/** Create a new node derived from this one by one of the eight
	* symmetry operations of the square.
	* @param mode a value from 0 to 7 indicative of which symmetry operation
	*  of the square is used to place the imaged node.
	*/
	node congru(const int mode) const
	{
		/** parameters \c mode from 0 to 7 generate a new node that is
		* either the original one (mode=0), rotated around the (0,0) center
		* of the coordinates by multiples of 90 degrees (mode=1 to 4), or
		* created in a two-step process of first reflecting the node across the
		* main (x=y) diagonal of the coordinates and then rotating it to 
		* a place away a multiple of 90 degrees.
		*/
		switch(mode)
		{
		case 0 :
			return *this ;
			break ;
		case 1 :
			// 90 degrees rotated
			return node(-xy[1],xy[0]) ;
			break ;
		case 2 :
			// 180 degrees rotated
			return node(-xy[0],-xy[1]) ;
			break ;
		case 3 :
			// 270 degrees rotated
			return node(xy[1],-xy[0]) ;
			break ;
		case 4 :
			// mirror main diagonal
			return node(xy[1],xy[0]) ;
			break ;
		case 5 :
			// mirror main diagonal then 90 degrees rotated
			return node(-xy[0],xy[1]) ;
			break ;
		case 6 :
			// mirror main diagonal then 180 degrees rotated
			return node(-xy[1],-xy[0]) ;
			break ;
		case 7 :
			// mirror main diagonal then 270 degrees rotated
			return node(xy[0],-xy[1]) ;
			break ;
		}
	}
protected:
private:
} ;

/** Node ordering operator.
* This introduces a ordering of a list of nodes by ordering with respect
* to the y coordinate (smaller y is interpreted as smaller node), and with
* respect to the x coordinate if the y coordinates are equal.
* @param n1 the node left from the minor sign.
* @param n2 the node right from the minor sign.
*/
bool operator< (const node &n1, const node &n2)
{
	if (n1.xy[1] < n2.xy[1])
		return true ;
	else if (n1.xy[1] > n2.xy[1])
		return false ;
	else if (n1.xy[0] < n2.xy[0])
		return true ;
	else
		return false ;
}


/** Node ordering operator.
* Two nodes are equal if they have the same x and the same y coordinate.
* @param n1 the node left from the equal sign.
* @param n2 the node right from the equal sign.
*/
bool operator== (const node &n1, const node &n2)
{
	return (n1.xy[0] == n2.xy[0] && n1.xy[1] == n2.xy[1]) ;
}

/** The edge class is a pair of nodes connected by an edge.
* The connectivity is implicit, and only the two coordinates of the
* pair of nodes are stored.
*/
class edge {
public:
	/* The pair of nodes, assumed  ordered such that nodes[0] < nodes[1].
	* To form an edge, nodes[0]!=nodes[1] in all cases.
	*/
	node nodes[2] ;

	/** The constructor with a pair of nodes.
	* @param n1 one of the two nodes.
	* @param n2 the other one of the two nodes.
	*/
	edge(const node n1, const node n2)
	{
		/* To ensure correct ordering, we call an interface function
		* of this class which swaps the two arguments if necessary 
		*/
		init(n1,n2) ;
	}

	/** A constructor with a node and an indicator of a direction
	* in which the other node is found.
	* @param n one node to be used for the edge.
	* @param dir a direction indicator from 0 to 3. This ressembles a
	*   wind rose indicator for the four potential directions on the
	*   square lattice.
	*/
	edge(const node n, int dir)
	{
#ifdef OLD
		/* This is in essence the slow version of the construction
		* since it involves an additional comparison with the init() function.
		*/
		switch(dir)
		{
		case 0:
			init(n, node(n.xy[0]+1,n.xy[1])) ;
			break ;
		case 1:
			init(n, node(n.xy[0],n.xy[1]+1)) ;
			break ;
		case 2:
			init(n, node(n.xy[0]-1,n.xy[1])) ;
			break ;
		case 3:
			init(n, node(n.xy[0],n.xy[1]-1)) ;
			break ;
		}
#else
		/* For dir values from 0 to 3 we create the pairing node in
		* the East, North, West and South direction (+x,+y,-x,-y) from the node \c n that
		* was provided explicitly.
		*/
		switch(dir)
		{
		case 0:
			nodes[0] = n ;
			nodes[1] = node(n.xy[0]+1,n.xy[1]) ;
			break ;
		case 1:
			nodes[0] = n ;
			nodes[1] = node(n.xy[0],n.xy[1]+1) ;
			break ;
		case 2:
			nodes[1] = n ;
			nodes[0] = node(n.xy[0]-1,n.xy[1]) ;
			break ;
		case 3:
			nodes[1] = n ;
			nodes[0] = node(n.xy[0],n.xy[1]-1) ;
			break ;
		}
#endif
	}

	/** This constructs a congruent version of the existing node.
	* @param mode a value form 0 to 7 used in nodes.congru() representing
	*   any of the eight ways (including the indentity/copy) to generate an
	*   image of the node by rotation and reflection.
	* @return the new edge rotated and/or reflected with respect to the
	*   origin.
	*/
	edge congru(const int mode) const
	{
		return edge( nodes[0].congru(mode), nodes[1].congru(mode)) ;
	}

	/** Translate the edge.
	* On return, the edge is moved by \c dx and dy into the two directions on the square grid.
	* @param dx the integer value for motion into the x direction
	* @param dy the integer value for motion into the y direction
	*/
	void shift(short dx, short dy)
	{
		nodes[0].xy[0] += dx ;
		nodes[1].xy[0] += dx ;
		nodes[0].xy[1] += dy ;
		nodes[1].xy[1] += dy ;
	}
protected:
private:
	/** Copy throug the nodes ensuring that the smaller node is placed first.
	* @param n1 the first of the nodes to be used for the new edge
	* @param n2 the second of the nodes to be used for the new edge
	*/
	void init(const node & n1, const node & n2)
	{
		/* The node comparison operator is involved to ensure a consistent
		* ordering of the two nodes in the short vector of nodes within the edge.
		*/
		if ( n1 < n2)
		{
			nodes[0]=n1;
			nodes[1]=n2;
		}
		else
		{
			nodes[0]=n2;
			nodes[1]=n1;
		}
	}
} ;

/** Edge ordering operator.
* The edge with the smaller first node is considered the smaller of the
* two edges. If the first node is shared between the two edges, the
* comparison is done on the second node.
* @param e1 the first of the two edges to be compared.
* @param e2 the second of the two edges to be compared.
* @return true if the first edge is considered the smaller one.
*/
bool operator< (const edge &e1, const edge &e2)
{
	if ( e1.nodes[0] < e2.nodes[0])
		return true ;
	else if ( e1.nodes[0] == e2.nodes[0])
	{
		if ( e1.nodes[1] < e2.nodes[1])
			return true ;
		else
			return false ;
	}
	else
		return false ;
}

/** Edge ordering operator.
* Edges are the same if the nodes are the same.
* @param e1 the first of the two edges to be compared.
* @param e2 the second of the two edges to be compared.
* @return true if the edges are the same.
*/
bool operator== (const edge &e1, const edge &e2)
{
	/* Since edges are sorted in their nodes, it's sufficient to compare the first and the
	* second nodes separately, not all four pairs of them... This is an attempt to speed
	* up the program.
	*/
	return ( e1.nodes[0] == e2.nodes[0] && e1.nodes[1] == e2.nodes[1]) ;
}

/** A tree is a list/vector of edges, all of which are different, and which
* share nodes such that one can transverse the graph on a trail that
* may reach each of the nodes..
*/
class tree {
public:
	/* The vector of all edges is the only full representation of the tree. */
	vector<edge> edges ;
	/* The \c bbox contains a bounding box, such that bbox[0] is the smallest x,
	* bbox[1] the smallest y, bbox[2] the largest x and bbox[3] the largest y value of
	* all nodes of the \c edges.
	*/
	short bbox[4] ;

	/* The default constructor creates an empty tree.
	* It is useful to generated vectors of trees and has no other significance.
	*/
	tree()
	{
	}

	/** The function congru() created a congruental image of the current tree
	* @param mode a number from 0 to 7 in the sense of edge::congru() which
	*   specifies which of the 8 symmetry operations (mode=0 is a copy) of the square
	*   is used to generate the new tree.
	* @return the rotated/reflected tree
	*/
	tree congru(const int mode) const
	{
		tree resul ;
		/* The algorithm is simply to create the image of each edge separately,
		* and add it to the new tree with addedge().
		*/
		for(int e=0; e < int(edges.size()); e++)
		{
			const edge cedge = edges[e].congru(mode) ;
			resul.addedge(cedge) ;
		}
		/* The typical rotation would have created an image in any quadrant of
		* the coordinate system. For quicker comparison of congruent pairs (which are
		* to be eliminated according to the algorithm), we translate each tree such
		* it is placed in the upper-right quadrant with the smallest x and the smallest
		* y coordinate of all nodes being both 0.
		*/
		resul.normalize() ;
		return resul ;
	}

	/** Create one more edge in the tree. This is done without
	* any further check such as ensuring that the edge is connected to the
	* already existing truncs or such as ensuring it is not yet part of the tree.
	* @param e the edge to be added
	*/
	void addedge(const edge & e)
	{
		edges.push_back(e) ;
	}

	/** Investigate whether the two trees are congruent.
	* The result is true if all the edges of this instance and the other instance
	* can be placed on each other by translation, rotation or mirror operation.
	* @param oth the other tree that is to be compared with the current instance.
	*   It is assumed that its \c bbox parameter is up-to-date.
	* @return true if this tree and the other one are congruent.
	*/
	bool iscongru(const tree & oth)
	{
		/* To avoid some overhead during the comparison, we compared some
		* rough parameters first: if the numbers of edges in the two trees differ,
		* the trees differ.
		*/
		if ( edges.size() != oth.edges.size() )
			return false ;
		/* The congruent operations of rotations and reflections can swap
		* the bounding box x and y coordinates, but not change them otherwise.
		* So we conclude that the two trees differ if this test of matching fails.
		*/
		if ( bbox[2] != oth.bbox[2] && bbox[2] != oth.bbox[3])
			return false ;
		if ( bbox[3] != oth.bbox[2] && bbox[3] != oth.bbox[3])
			return false ;

		/* The final more constly operation of deciding on congruence is to create
		* all 8 congruent images of \c oth and compare each of them with the
		* current instance. To speed up the comparison, we sort the edge list
		* of the current instance and each edge list of the 8 images created of the other
		* instance. This ought result in a comparison time which is only linear in the
		* number of edges.
		*/
		sort(edges.begin(),edges.end() ) ;

		/* This is the loop over all possible 8 congruental versions of the other tree. */
		for(int c=0; c < 8 ;c++)
		{
			/* The temporary mage copy of the other tree by calling congru(). */
			tree imag = oth.congru(c) ; 
			/* Sorting the edge list of the image copy and figuring out whether the
			* image lists are then the same is done with two calls of the STL library.
			* If we found that any of the 8 versions equals the current instance,
			* we do not need to investigate the other versions and may conclude
			* that both trees are congruental.
			*/
			sort(imag.edges.begin(),imag.edges.end() ) ;
			if ( equal(edges.begin(),edges.end(),imag.edges.begin()) )
				return true ;
		}
		/* If we have fallen through the loop over the 8 images, none of these
		* was congruental and the result is a 'no'.
		*/
		return false ;
	}

	/** A function to generate all trees with one more edge than the existing instance.
	* The current instance is the parent from which all possible ways of adding
	* one more edge to any existing node are tested. The children are reduced to
	* a set of mutually incongruental trees and returned in a vector.
	* @return a vector of the incongruental trees generated.
	*/
	vector<tree> childs()
	{
		/* \c resul collects all new trees with one more edge than the current tree */
		vector<tree> resul ;

		/* For some efficiency, we generate a list of all nodes of the existing tree,
		* which are typically shared by between 1 and 4 edges. The idea is to spawn
		* new edges starting from any of these into any of the four directions
		* admitted by the square lattice.
		*/
		const vector<node> thisnod = allnodes() ;
		//cout << thisnod.size() << " nodes " << edges.size() << " edges\n" ;

		/* In a loop over all existing nodes we can generate all child trees
		* that have the same skeleton as the current one, but one more edge.
		* This works since a tree is a connected graph, so edges cannot be isolated.
		*/
		for(int n=0; n < int(thisnod.size()) ; n ++)
		{
			/* For each of the four directions of the node we test whether this
			* can be the directino of a new edge for a potentially new child tree.
			*/
			for(int dir=0; dir <4 ; dir++)
			{
#ifdef OLD
				/* This is a more costly (slower) version of the algorithm below */
				const edge tste( thisnod[n],dir) ;
				if ( tstedge(tste) )
				{
					tree candi(*this) ;
					candi.addedge(tste) ;
					candi.normalize() ;
					bool congru=false ;
					for(int pres=0; pres < int(resul.size()); pres++)
					{
						if ( resul[pres].iscongru(candi) )
						{
							congru= true ;
							break ;
						}
					}
					if ( !congru)
						resul.push_back(candi) ;
				}
#else
				/* For this particular node and direction, we test whether
				* the node and the new node are a candidate for extension. We generate
				* a new node into the particular direction, and admit this as a new
				* edge if this new site is not yet already occupied.
				*/
				const node tste( thisnod[n],dir) ;
				/* Test on occupation of the point in the square lattice uses havnode(). */
				if ( !havnode(tste) )
				{
					/* If this node is not yet occupied (allowed since a tree is cycle-free),
					* we create the new child by first copying this tree, and adding
					* a new edge that spreads between the current node and the new one.
					*/
					tree candi(*this) ;
					candi.addedge(edge(thisnod[n],tste)) ;

					/* Adding the edge may have extended the new tree into the
					* quadrant to the left or below, so we shift the tree to the right
					* or up with normalize() to place it into the first quadrant
					* of the coordinate system.
					*/
					candi.normalize() ;
					/* This new child tree may already exist in the list of trees
					* generated before. We test wheter it is congruent to any fo the
					* child trees already accumulated so far in the result vector.
					*/
					bool congru=false ;
					for(int pres=0; pres < int(resul.size()); pres++)
					{
						if ( resul[pres].iscongru(candi) )
						{
							congru= true ;
							break ;
						}
					}
					/* If no version of the new child tree is found in the
					* current list of children, we put it into the result list.
					*/
					if ( !congru)
						resul.push_back(candi) ;
				}
#endif
			}
		}

		/* The result list of children that are mutualy not congruental is returned */
		return resul ;
	}
protected:
private:

	/* The bounding box \c bbox is updated to reflect the actual
	* size of the tree.
	*/
	void updateBB()
	{
		/* The algorithm scans each edge and both nodes in the edge
		* for the maximum and minimum coordinate in x and y, and updates
		* the four bounding box parameters where the edge sprawls outside
		* the bounding box seen so fare.
		*/
		if ( edges.size() )
		{
			bbox[2] = bbox[0] = edges[0].nodes[0].xy[0] ;
			bbox[3] = bbox[1] = edges[0].nodes[0].xy[1] ;
		}
		for(int e=0; e < int(edges.size()); e++)
		{
			bbox[0] = min(bbox[0],edges[e].nodes[0].xy[0]) ;
			bbox[0] = min(bbox[0],edges[e].nodes[1].xy[0]) ;
			bbox[1] = min(bbox[1],edges[e].nodes[0].xy[1]) ;
			// bbox[1] = min(bbox[1],edges[e].nodes[1].xy[1]) ; nodes are sorted: sufficient to look at nodes[0] for the minimum y-component
			bbox[2] = max(bbox[2],edges[e].nodes[0].xy[0]) ;
			bbox[2] = max(bbox[2],edges[e].nodes[1].xy[0]) ;
			//bbox[3] = max(bbox[3],edges[e].nodes[0].xy[1]) ; nodes are sorted: sufficient to look at nodes[1] for the maximum y
			bbox[3] = max(bbox[3],edges[e].nodes[1].xy[1]) ;
		}
	}

	/** Normalization means rigid translation of the tree coordinates
	* (of all nodes in all edges) such that the tree is tightly fitting into
	* the first quadrant and the minimum x and y coordinates over all nodes are zero and zero.
	*/
	void normalize()
	{
		updateBB() ;
		if ( bbox[0] || bbox[1] )
		{
			const short dx = -bbox[0] ;
			const short dy = -bbox[1] ;
			for(int e=0; e < int(edges.size()); e++)
				edges[e].shift(dx,dy) ;
			bbox[0] =bbox[1] = 0 ;
			bbox[2] += dx ;
			bbox[3] += dy ;
		}
	}

	/* Test whether \c n is already part of the tree.
	* @param n the node to be found or not found in the current instance.
	* @return true if the node is part of any edge in the tree.
	*/
	bool havnode(const node & n) const
	{
		/* In a simple linear search strategy, all edges and both nodes of each
		* of them is compared with the edge \c n.
		*/
		for(int e=0; e < int(edges.size()); e++)
			if ( edges[e].nodes[0] == n || edges[e].nodes[1] == n)
				/* If we found a node in the tree that matches the candidate,
				* we report this as early as possible for efficiency.
				*/
				return true ;

		/* At this point of the program, we have fallen through the double loop over
		* edges and nodes and have not found any match. So the result is a 'no'.
		*/
		return false ;
	}

#ifdef OLD
	/** Test wheter 'e' is an acceptable/compatible choice for a new edge
	*/
	bool tstedge(const edge & e) const
	{
		const bool haven1 = havnode( e.nodes[0]) ;
		const bool haven2 = havnode( e.nodes[1]) ;
		/* return true if have node 0 but not node 1 or vice versa */
		return ( haven1 != haven2) ;
	}
#endif
	/** Test an edge as a candidate for a new edge in the current tree.
	* @param n one node of current instance of a tree.
	* @param dir a lookup direction (0 to 3) for extension along a new edge.
	* @return true (equalent to admitted) indicating that the new site
	*   for the partner node of \c n is not yet occupied by another node
	*   in the current tree.
	*/
	bool tstedge(const node & n, int dir) const
	{
		return !havnode(node(n,dir)) ;
	}

	/** Generate a list of all nodes in the current tree.
	* @return the vector of all nodes (duplicates removed) of the current tree.
	*/
	vector<node> allnodes() const
	{
		/** This is the result vector, empty at creation */
		vector<node> resul ;
		/* We linearly search through all edges that host two nodes each */
		for(int e=0; e < int(edges.size()) ;e++)
		{
			/* The algorithm::find() operation of the STL is used to compare
			* the node of the current each with the set of nodes already seen
			* by visiting the previous edges.
			*/
			vector<node>::iterator fi=resul.begin() ;
			vector<node>::iterator la=resul.end() ;
			if ( find(fi,la,edges[e].nodes[0]) == resul.end() )
				resul.push_back( edges[e].nodes[0]) ;
			fi=resul.begin() ;
			la=resul.end() ;
			if ( find(fi,la,edges[e].nodes[1]) == resul.end() )
				resul.push_back( edges[e].nodes[1]) ;
		}
		return resul ;
	}
} ;

/** Tree comparison operator.
* @param t1 the first tree, supposed to be normalized in the first quadrant.
* @param t2 the second tree, supposed to be normalized in the first quadrant
* @return true if the two trees overlap completely (share the same list of edges).
*/
bool operator== (tree &t1, tree &t2)
{
	/* Trees of differnt edge counts are different */
	if ( t1.edges.size() != t2.edges.size() )
		return false ;
	/* Trees of different bounding boxes are also different */
	else if ( t1.bbox[2] != t2.bbox[2] && t1.bbox[2] != t2.bbox[3]
			|| t1.bbox[3] != t2.bbox[2] && t1.bbox[3] != t2.bbox[3] )
		return false ;
	else
	{
		/* Trees of same edge count and bounding box are the same if the lists
		* of edges, after sorting according to the same criteria of comparison operators,
		* are the same.
		*/
		sort(t1.edges.begin(),t1.edges.end() ) ;
		sort(t2.edges.begin(),t2.edges.end() ) ;
		return equal(t1.edges.begin(),t1.edges.end(),t2.edges.begin()) ;
	}
}

/** A human readable ASCII output version of the tree.
* A matrix of dots (empty places), oh's (nodes) and dashes (up or sideways)
* is plotted in ASCII art fashion.
* @param t the tree to be printed.
* @param os the output stream to be used.
*/
ostream & operator<<(ostream & os, const tree &t)
{
	const int xdim = t.bbox[2]+1 ;	// bbox only indicates the hightes occupied member
	const int ydim = t.bbox[3]+1 ;
	char **arr = new char *[2*ydim-1] ;
	for(int r=0 ; r < 2*ydim-1; r++)
	{
		arr[r] = new char [2*xdim-1] ;
		for(int c=0; c < 2*xdim-1; c++)
			arr[r][c]='.' ;
	}
	for(int e=0; e < int(t.edges.size()) ; e++)
	{
		int x0= t.edges[e].nodes[0].xy[0] ;
		int y0= t.edges[e].nodes[0].xy[1] ;
		arr[2*y0][2*x0]='o' ;
		int x1= t.edges[e].nodes[1].xy[0] ;
		int y1= t.edges[e].nodes[1].xy[1] ;
		arr[2*y1][2*x1]='o' ;
		if ( x0 == x1)
			arr[y0+y1][2*x0]='|' ;
		else
			arr[2*y0][x0+x1]='-' ;
	}
	for(int r=0; r < 2*ydim-1 ; r++)
	{
		for(int c=0; c < 2*xdim-1; c++)
			os << arr[r][c] ;
		os << endl ;
	}
	os << endl ;

	for(int r=0 ; r < 2*ydim-1; r++)
		delete [] arr[r] ;
	delete [] arr ;
	return os ;
}

/** Create a new generation of child trees by generation of all in-congruental
* childs of the trees in \c paren.
* @param paren the existing forest of parent trees
* @param verbose if true the new child trees are displayed on stdout.
* @return the vector of all incongruental trees with one more edge than
*  those of the \c paren list.
*/
vector<tree> gene(vector<tree> paren, bool verbose)
{
	/* The result vector */
	vector<tree> resul ;
	/* New trees are generated by spawning from each of the parent trees in turn, in a loop */
	for(int t=0; t < int(paren.size()) ; t++)
	{
		/* The child trees of this particular parent are collected in an intermediate
		* forest \c chi.
		*/
		vector<tree> chi = paren[t].childs() ;
		for(int c=0 ; c < chi.size() ; c++)
		{
#if 0
			cout << "--------------\n" ;
			cout << paren[t] << endl ;
			cout << "generated\n" ;
			cout << chi[c] << endl ;
			cout << "--------------\n" ;
#endif
			/* Each of this child trees is tested for congruency relative to the
			* existing full set of trees in the result vector.
			*/
			bool congru=false ;
			for(int pres=0; pres < int(resul.size()); pres++)
			{
				if ( resul[pres].iscongru(chi[c]) )
				{
					congru= true ;
					break ;
				}
			}
			/* If a congruent version is not yet there, we add the new child tree,
			* and if the verbose option was used, show and enumerate it on stdout.
			*/
			if ( !congru)
			{
				if( verbose )
				{
					cout << resul.size() << ":\n" ;
					cout << chi[c] << endl ;
				}
				resul.push_back(chi[c]) ;
			}
		}
	}
	return resul ;
}

/** The main program.
* There is only one command line option, -v, which means that
* an explicit list of each generatino of trees is plotted in ASCII art
* on stdout. Otherwise, only the raw counts are reported.
*/
int main(int argc, char *argv[])
{

	/* We test for the presence of the -v option */
	bool verbose =false ;
	char oc ;
	while (  (oc=getopt(argc,argv,"v")) != -1 )
	{
		switch(oc)
		{
		case 'v' :
			verbose = true ;
			break ;
		case '?' :
			cerr << "Invalid command line option " << optopt << endl ;
			break ;
		}
	}

	/* The fundamental root tree (great grand root tree) is the one from (0,0) to (1,0),
	* with one edge and two nodes only.
	*/
	tree root ;
	edge strt(node(0,0),node(1,0)) ;
	root.addedge(strt) ;
	root.bbox[0]=0 ;root.bbox[1]=0 ;
	root.bbox[2]=1 ;root.bbox[3]=0 ;

	vector<tree> lvl1 = root.childs() ;
#if 0
	for(int t=0 ;t < int(lvl1.size()) ; t++)
	{
		cout << t << ":\n" ;
		for(int e=0; e < lvl1[t].edges.size() ; e++)
		{
			const edge & ee = lvl1[t].edges[e] ;
			cout << "from " << ee.nodes[0].xy[0] << " " << ee.nodes[0].xy[1] << " to " ;
			cout << ee.nodes[1].xy[0] << " " << ee.nodes[1].xy[1] << endl ;
		}
	}
#endif
	cout << "total " << lvl1.size() << endl ;

	/* The program generates the first level generation from all
	* children of the root tree, the second level generatio from all children of
	* all trees in the first level etc. Obviously one can extend this further
	* beyond the hard-coded maximum of 8 generations of trees.
	*/
	vector<tree> lvl2 = gene(lvl1,verbose) ;
	cout << "total " << lvl2.size() << endl ;

	vector<tree> lvl3 = gene(lvl2,verbose) ;
	cout << "total " << lvl3.size() << endl ;

	vector<tree> lvl4 = gene(lvl3,verbose) ;
	cout << "total " << lvl4.size() << endl ;

	vector<tree> lvl5 = gene(lvl4,verbose) ;
	cout << "total " << lvl5.size() << endl ;

	vector<tree> lvl6 = gene(lvl5,verbose) ;
	cout << "total " << lvl6.size() << endl ;

	vector<tree> lvl7 = gene(lvl6,verbose) ;
	cout << "total " << lvl7.size() << endl ;

	vector<tree> lvl8 = gene(lvl7,verbose) ;
	cout << "total " << lvl8.size() << endl ;

	return 0 ;
}