/*************************************************************** 
	The GPDD Project (GPDD1) 
 
	GPDD stands for "Graph-Pair Decision Diagram" 
 
Contributors: 
 
   CHEN, Weiwei			(MS 2007 graduate) 
   MENG, Xiaoxuan		(MS 2009 graduate) 
   ZENG, Mei			(MS 2010 graduate) 
   LI, Ji				(MS 2010 graduate) 
   MA, Diming			(MS 2011 graduate) 
 
   SHI, Guoyong 		(Prof)  shiguoyong@ic.sjtu.edu.cn 
 
   School of Microelectronics 
   Shanghai Jiao Tong University 
   (c) 2005 - 2010 
******************************************************************/ 



/**************************************************************
  GPDD.cpp
  GPDD is the main class for Graph-Pair Decision Diagram.

  The two kernel functions for graph reduction are
	GPDD::CreateLeftGraph() and GPDD::CreateRightGraph().


	In this implementation, the following components have been removed because they
	are not very important for the purpose of the Sensitivity package. 
	By removing these insignificant functions, the GPDD package is greatly condensed.
	(G. Shi, 08-31-2010)

	(1) The approximate H(s) analysis functionalities.
	(2) The redundancy reduction during GPDD construction functionalities (those function
	     names carrying an 'S', like ReduceS(), GCS(), and doGCS(), etc.
	(3) The function of GPDD construction with redundancy removal on-the-fly,
	     which uses GCS() and ReduceS(), is also removed.

***************************************************************/

#include <cstdlib>
#include <cstring>
#include <map>
#include <iostream>
#include <iomanip>
#include <stack>
#include <string>

#include "gpdd_defs.h"
#include "gpdd_global.h"
#include "utility.h"
#include "GPDD.h"

using namespace std;

//extern int compare (const void*, const void*);
//extern int compare_g (const void*, const void*);
//extern int compare_l (const void*, const void*);

extern const char *EdgeTypeName[];

GPDD::GPDD() :
 	AnyObject(), CountOfMinusTerms(0), CountOfPlusTerms(0), TotalExpandedGPDDnodes(0),
	NodeCreationCount(0), SubgraphCount(0), SubgraphSharingCount(0),
 	GPDDnodeSharingCount(0), TotalGPDDnodes(0), NumNodesRemoved(0),
	PermMark(0), NodeMark(0), MaxNumNodes(0)
{
	cout << LINE_INFO << "In constructor -" << endl;
	
	NumSymbols = 0;
	NumTerms = 0;
	
	EdgeList				= NULL;
	LumpedSymbList		= NULL;
	ExpSymbList			= NULL;

	NumEdges	= 0;
	NumEdges_L	= 0;
	NumEdges_R	= 0;

	NumNodes			= 0;
	VsN					= 0;
	VcN					= 0;
	CsN					= 0;
	CcN					= 0;
	DpenSrcN			= 0;
	StackTop				= 0;

	pGPDDroot			= NULL;
	VertexOne			= NULL;
	VertexZero			= NULL;
	PathStack			= NULL;
	NumNodesRemaining		= 0;
	sCoeffNodeCount	= 0;

	sCoeffOne			= NULL;
	sCoeffZero			= NULL;

	ExpGPDDroot			= NULL;
	ExpOne			= NULL;
	ExpZero			= NULL;

	LumpedExpandedFlag = false;
	sExpandedFlag = false;
	TraversalMark = false;
	
}


#if	00
/**
	Constructor for HybridSim.
*/
GPDD::GPDD(int nn, int en, vector<Device *> &dev_vec) :
 	AnyObject(), TotalNodes_(0), CountOfMinusTerms(0), CountOfPlusTerms(0), TotalExpandedGPDDnodes(0),
	NodeCreationCount(0), SubgraphCount(0), SubgraphSharingCount(0),
	GPDDnodeSharingCount(0), TotalGPDDnodes(0), NumNodesRemoved(0),
	PermMark(0), NodeMark(0), MaxNumNodes(0)
{
	cout << LINE_INFO << "Construct by the device vector " << endl;
	
	int			i, j;

	NumNodes		= nn;
	NumEdges		= en;
	NumSymbols			= 0;
	NumTerms			= 0;
	EdgeList			= NULL;
	LumpedSymbList			= NULL;
	SymbListTemp	= NULL;
	ExpSymbList		= NULL;
	VsN				= 0;
	VcN				= 0;
	CsN				= 0;
	CcN				= 0;
	DpenSrcN		= 0;
	StackTop			= 0;
	
	SymbOne.symbindex = -1;
	SymbOne.symbname = CopyStr("One");
	
	SymbZero.symbindex = -2;
	SymbZero.symbname = CopyStr("Zero");
	
	VertexOne = new GPDDnode(&SymbOne, +1, +1, NULL, NULL);
	VertexZero = new GPDDnode(&SymbZero, +1, +1, NULL, NULL);
	
	PathStack		= NULL;
	
	ExpGPDDroot			= NULL;

	uniqSubgraph			= new HashTree*[SUBGRAPH_PRIME];
	uniqGPDD				= new HashTree**[GPDD_PRIME_LEVEL];
	uniqCoeffGPDD			= new HashTree**[TAG_PRIME1];
	uniqExpGPDD			= new HashTree**[TAG_PRIME1];


	for(i = 0; i < SUBGRAPH_PRIME; i ++)
	{
		uniqSubgraph[i] = NULL;
	}

	for(i = 0; i < TAG_PRIME1; i++)
	{
		uniqExpGPDD[i] = new HashTree*[TAG_PRIME2];
		for(j = 0; j < TAG_PRIME2; j ++)
			uniqExpGPDD[i][j] = NULL;
	}/**/

	for(i = 0; i < GPDD_PRIME_LEVEL; i ++)
	{
		uniqGPDD[i] = new HashTree*[GPDD_PRIME];
		for(j = 0; j < GPDD_PRIME; j ++)
			uniqGPDD[i][j] = NULL;
	}
	
#ifdef HASHSCOEFFNODE 
	for(i = 0; i < TAG_PRIME1; i ++)
	{
		uniqCoeffGPDD[i] = new HashTree*[TAG_PRIME2];
		for(j = 0; j < TAG_PRIME2; j ++)
			uniqCoeffGPDD[i][j] = NULL;
	}
#endif

	if ( !CreateEdgeList(dev_vec) )
	{
		cerr << "\nEnumerator(): Failed to create Edge List." 
	  		 << "\nProgram terminated." << endl;
		exit(1);
	}
	if ( CreateSymbList() == STOP)
	{
		cerr << "\nEnumerator(): Failed to create Symbol List." 
	  		 << "\nProgram terminated." << endl;
		exit(1);
	}

  	srand((int) time((time_t*) NULL));
 	// My_ZRand = rand();

	PermMark = new int[NumNodes];
	NodeMark = new int[NumNodes];
}
#endif


/**
	Constructor without special treatment of MOSFET.
*/
GPDD::GPDD(int nn, int en, Edge_g *edge_list) :
 	AnyObject(), CountOfMinusTerms(0), CountOfPlusTerms(0), TotalExpandedGPDDnodes(0),
	NodeCreationCount(0), SubgraphCount(0), SubgraphSharingCount(0),
	GPDDnodeSharingCount(0), TotalGPDDnodes(0), NumNodesRemoved(0),
	PermMark(0), NodeMark(0), MaxNumNodes(0)
{
	cout << LINE_INFO << "Construction by the edge list " << endl;
	
	int		i, j;

	NumNodes = nn;
	NumEdges = en;
	NumEdges_L	= 0;
	NumEdges_R	= 0;
	
	NumSymbols = 0;
	NumTerms = 0;
	EdgeList = NULL;
	LumpedSymbList = NULL;
	SymbListTemp	= NULL;
	ExpSymbList = NULL;
	
	VsN = 0;
	VcN = 0;
	CsN = 0;
	CcN = 0;
	DpenSrcN = 0;
	StackTop	 = 0;

	
	SymbOne.symbindex = -1;
	SymbOne.symbname = CopyStr("One");
	
	SymbZero.symbindex = -2;
	SymbZero.symbname = CopyStr("Zero");
	
	VertexOne = new GPDDnode(&SymbOne, +1, +1, NULL, NULL);
	VertexZero = new GPDDnode(&SymbZero, +1, +1, NULL, NULL);
	
	PathStack = NULL;
	
	ExpGPDDroot = NULL;

	uniqSubgraph = new HashTree*[SUBGRAPH_PRIME];
	uniqGPDD	= new HashTree**[GPDD_PRIME_LEVEL];
	uniqCoeffGPDD = new HashTree**[TAG_PRIME1];
	uniqExpGPDD = new HashTree**[TAG_PRIME1];


	for(i = 0; i < SUBGRAPH_PRIME; i ++)
	{
		uniqSubgraph[i] = NULL;
	}

	for(i = 0; i < TAG_PRIME1; i++)
	{
		uniqExpGPDD[i] = new HashTree*[TAG_PRIME2];
		for(j = 0; j < TAG_PRIME2; j ++)
			uniqExpGPDD[i][j] = NULL;
	}

	for(i = 0; i < GPDD_PRIME_LEVEL; i ++)
	{
		uniqGPDD[i] = new HashTree*[GPDD_PRIME];
		for(j = 0; j < GPDD_PRIME; j ++)
			uniqGPDD[i][j] = NULL;
	}
	
#ifdef HASH_S_COEFF_NODE 
	for(i = 0; i < TAG_PRIME1; i ++)
	{
		uniqCoeffGPDD[i] = new HashTree*[TAG_PRIME2];
		for(j = 0; j < TAG_PRIME2; j ++)
			uniqCoeffGPDD[i][j] = NULL;
	}
#endif

	if ( !CreateEdgeList(edge_list) )
	{
		cerr << LINE_INFO << "\nEnumerator(): Failed to create Edge List." 
	  		 << "\nProgram terminated." << endl;
		exit(1);
	}
	if ( CreateSymbList() == STOP )
	{
		cerr << LINE_INFO << "\nEnumerator(): Failed to create Symbol List." 
	  		 << "\nProgram terminated." << endl;
		exit(1);
	}

  	srand((int) time((time_t*) NULL));
 	// My_ZRand = rand();

	PermMark = new int[NumNodes];
	NodeMark = new int[NumNodes];

	LumpedExpandedFlag = false;
	sExpandedFlag = false;
	TraversalMark = false;

}

/**
	Constructor with MOSFET.
*/
GPDD::GPDD(int num_nodes,  int num_edges, int num_mosfets, Edge_g *edge_list) :
 	AnyObject(), CountOfMinusTerms(0), CountOfPlusTerms(0), TotalExpandedGPDDnodes(0),
	NodeCreationCount(0), SubgraphCount(0), SubgraphSharingCount(0),
	GPDDnodeSharingCount(0), TotalGPDDnodes(0), NumNodesRemoved(0),
	PermMark(0), NodeMark(0), MaxNumNodes(0)
{
	cout << LINE_INFO << "Construction with MOSFET" << endl;
	
	NumMOSDevices = num_mosfets;
//	MOS_List = NULL;
	
	int		i, j;

	NumNodes = num_nodes;
	NumEdges = num_edges;
	NumEdges_L	= 0;
	NumEdges_R	= 0;

	NumSymbols = 0;
	NumTerms = 0;
	EdgeList = NULL;
	LumpedSymbList = NULL;
	SymbListTemp	= NULL;
	ExpSymbList = NULL;
	
	VsN = 0;
	VcN = 0;
	CsN = 0;
	CcN = 0;
	DpenSrcN = 0;
	StackTop	 = 0;

	
	SymbOne.symbindex = -1;
	SymbOne.symbname = CopyStr("One");
	
	SymbZero.symbindex = -2;
	SymbZero.symbname = CopyStr("Zero");
	
	VertexOne = new GPDDnode(&SymbOne, +1, +1, NULL, NULL);
	VertexZero = new GPDDnode(&SymbZero, +1, +1, NULL, NULL);
	
	PathStack = NULL;
	
	ExpGPDDroot = NULL;

	uniqSubgraph = new HashTree*[SUBGRAPH_PRIME];
	uniqGPDD	= new HashTree**[GPDD_PRIME_LEVEL];
	uniqCoeffGPDD = new HashTree**[TAG_PRIME1];
	uniqExpGPDD = new HashTree**[TAG_PRIME1];


	for(i = 0; i < SUBGRAPH_PRIME; i ++)
	{
		uniqSubgraph[i] = NULL;
	}

	for(i = 0; i < TAG_PRIME1; i++)
	{
		uniqExpGPDD[i] = new HashTree*[TAG_PRIME2];
		for(j = 0; j < TAG_PRIME2; j ++)
			uniqExpGPDD[i][j] = NULL;
	}

	for(i = 0; i < GPDD_PRIME_LEVEL; i ++)
	{
		uniqGPDD[i] = new HashTree*[GPDD_PRIME];
		for(j = 0; j < GPDD_PRIME; j ++)
			uniqGPDD[i][j] = NULL;
	}
	
#ifdef HASH_S_COEFF_NODE 
	for(i = 0; i < TAG_PRIME1; i ++)
	{
		uniqCoeffGPDD[i] = new HashTree*[TAG_PRIME2];
		for(j = 0; j < TAG_PRIME2; j ++)
			uniqCoeffGPDD[i][j] = NULL;
	}
#endif

	if ( CreateEdgeList(edge_list) != OK )
	{
		cerr << LINE_INFO << "\n\tFailed to create Edge List. Program is terminated." 
	  		<< endl;
		exit(1);
	}

	cout << LINE_INFO << "Run to here." << endl << flush;

	if ( CreateSymbList() == STOP )
	{
		cerr << LINE_INFO << "\n\tFailed to create Symbol List. Program is terminated." 
	  		<< endl;
		exit(1);
	}

  	srand((int) time((time_t*) NULL));
 	// My_ZRand = rand();

	PermMark = new int[NumNodes];
	NodeMark = new int[NumNodes];

	LumpedExpandedFlag = false;
	sExpandedFlag = false;
	TraversalMark = false;

}



GPDD::~GPDD()
{
	int i, j;

	if (PermMark)
	{
		delete[] PermMark;
		PermMark = 0;
	}
	if (NodeMark)
	{
		delete[] NodeMark;
		NodeMark = 0;
	}

#ifdef CALCULATE_S_COEFFICIENTS
	CleanExpandedGPDDnodes();
#endif

	if (EdgeList)
	{
		delete[] EdgeList;
		EdgeList = 0;
	}

	if (VertexOne)
	{
		delete VertexOne;
		VertexOne = 0;
	}

	if (VertexZero)
	{
		delete VertexZero;
		VertexZero = 0;
	}

	if (PathStack)
	{
		delete[] PathStack;
		PathStack = 0;
	}

	if (LumpedSymbList)
	{
		SymbNode	*pNode1, *pNode2;
		for (i = 0; i < NumSymbols; i++)
		{
			pNode2 = LumpedSymbList[i]->nextP;
			while (pNode2)
			{
				pNode1 = pNode2;
				pNode2 = pNode2->next;
				delete[] pNode1->symbname;
				delete pNode1;
			}
			delete[] LumpedSymbList[i]->symbname;
			delete LumpedSymbList[i];
			LumpedSymbList[i] = 0;
		}
		delete [] LumpedSymbList;
		LumpedSymbList = 0;
	}
//	delete [] ExpSymbList;

	if (uniqSubgraph)
	{
		for(i = 0 ; i < SUBGRAPH_PRIME ; i ++ )
		{
			if(!uniqSubgraph[i])
				continue;
			uniqSubgraph[i]->SetDeleteFlag(TRUE,FALSE);    
			delete uniqSubgraph[i];
			uniqSubgraph[i] = 0;
		}
		delete []uniqSubgraph;
		uniqSubgraph = 0;
	}

	if (uniqGPDD)
	{
		for(i = 0 ; i < GPDD_PRIME_LEVEL ; i ++ )
		{
			if(!uniqGPDD[i])
				continue;
			for(j = 0; j < GPDD_PRIME; j ++)
			{
				if(!uniqGPDD[i][j])
					continue;
				uniqGPDD[i][j]->SetDeleteFlag(TRUE,FALSE);
				delete uniqGPDD[i][j];
				uniqGPDD[i][j] = 0;
			}
			delete []uniqGPDD[i];
			uniqGPDD[i] = 0;
		}
		delete []uniqGPDD;
		uniqGPDD = 0;
	}

	if (uniqExpGPDD)
	{
		for(i = 0; i < TAG_PRIME1; i++)
		{
			if(!uniqExpGPDD[i])
				continue;
			for(j = 0; j < TAG_PRIME2; j ++)
			{
				if(!uniqExpGPDD[i][j] )
					continue;
				uniqExpGPDD[i][j]->SetDeleteFlag(TRUE, FALSE);
				delete uniqExpGPDD[i][j];
				uniqExpGPDD[i][j] = 0;
			}
			delete []uniqExpGPDD[i];
			uniqExpGPDD[i] = 0;
		}
		delete []uniqExpGPDD;
		uniqExpGPDD = 0;
	}

	NodeList.clear();


#ifdef HASHSCOEFFNODE 
	for (i = 0; i < TAG_PRIME1; i ++)
	{
		if (uniqCoeffGPDD[i])
		{
			if (!uniqCoeffGPDD[i])
				continue;
			for (j = 0; j < TAG_PRIME2; j ++)
			{
				if (!uniqCoeffGPDD[i][j])
					continue;
				uniqCoeffGPDD[i][j]->SetDeleteFlag(TRUE, FALSE);
				delete uniqCoeffGPDD[i][j];
				uniqCoeffGPDD[i][j] = 0;
			}
			delete[] uniqCoeffGPDD[i];
			uniqCoeffGPDD[i] = 0;
		}
	}

	delete[] uniqCoeffGPDD;
	uniqCoeffGPDD = 0;
#endif

#ifdef CALCULATE_S_COEFFICIENTS
#if 0
	delete sCoeffOne;
	sCoeffOne = 0;
	delete sCoeffZero;
	sCoeffZero = 0;
#endif
	// sCoeffOne and sCoeffZero will be deleted
	// when deleting pSexpOne and ExpZero.
	delete ExpOne;
	ExpOne = 0;
	delete ExpZero;
	ExpZero = 0;
#endif

	cout << LINE_INFO << "Deconstructed" << endl;
}


#if 	00
/**
  Use the device list to create the list of edges (for the subcircuit graph).

  This function is for use in HIERARCHICAL SYMBOLIC ANALYSIS (the HYBRIDSIM package).
*/
int GPDD::CreateEdgeList(vector<Device *> &dev_vec)
{
	cout << LINE_INFO << "Using the list of devices." << endl;
	
    if (NumEdges <= 0)
    {
        cerr <<"\nNo edge in graph.";
        return 0;
    }
    else
    {
        EdgeList = new Edge[NumEdges];
    }

    int e  = 0;
    int max_n = 0;

	cout << LINE_INFO << "Listing the edges -- " << endl;

    vector<Device *>::iterator it = dev_vec.begin();
    for( ; it!=dev_vec.end(); it++)
    {
        EdgeList[e].node1 = (*it)->GetNodeIndex(0);
        EdgeList[e].node2 = (*it)->GetNodeIndex(1);
        EdgeList[e].lumped = 0;
        EdgeList[e].edgenum = e;
        EdgeList[e].value = (*it)->GetValue();

	//Sensitivity of the device to a certain parameter
	//EdgeList[e].Sens = (*it)->GetSens();

//#ifdef TRACE
	printf("\n    Edge(%d, %d): %s", EdgeList[e].node1, EdgeList[e].node2, (*it)->GetName().c_str());
//#endif

        // Check syntax.
        if (EdgeList[e].node1 == EdgeList[e].node2)
        {
            cout << "\nError in CreateEdgeList(): Edge #" << e 
                   << " (" << EdgeList[e].name << ") " << "is a self-loop.";
            delete [] EdgeList;
            EdgeList = 0;
            return 0;		// fail
        }
	cout << endl;

        string str = (*it)->GetName();

        char edge_name_prefix = str[0];
        edge_name_prefix = toupper(edge_name_prefix);

        //name of controlled source and controlling source
        string cName, sName;
        switch (edge_name_prefix)
        {
            case 'R':
                EdgeList[e].name = CopyStr((*it)->GetName().c_str());
                //EdgeList[e].value = (*it)->GetValue();
                EdgeList[e].type = Z;
                break;
            case 'C':
            case 'L': 
                EdgeList[e].name = CopyStr((*it)->GetName().c_str());
                //EdgeList[e].value = (*it)->GetValue();
                EdgeList[e].type = Y;
                break;
            case 'E':

                sName = (*it)->GetName();//"VS" + (*it)->GetName();
                cName = (*it)->GetName();//"VC" + (*it)->GetName();
		EdgeList[e].name = CopyStr(sName.c_str());
                EdgeList[e].pname = CopyStr(cName.c_str());
                //EdgeList[e].value = edge->value;
                EdgeList[e].type = VS;
                VsN ++;

                //Add VC
                e++;
                EdgeList[e].name = CopyStr(cName.c_str());
                EdgeList[e].pname = CopyStr(sName.c_str());

                EdgeList[e].node1 = (*it)->GetNodeIndex(2);
                EdgeList[e].node2 = (*it)->GetNodeIndex(3);
                EdgeList[e].lumped = 0;
                EdgeList[e].edgenum = e;
                EdgeList[e].type = VC;

                VcN++;

                EdgeList[e-1].pnum = e;
                EdgeList[e].pnum = e-1;
                break;
            case 'F':	
                sName = (*it)->GetName();//"CS" + (*it)->GetName();
                cName = (*it)->GetName();//"CC" + (*it)->GetName();
                EdgeList[e].name = CopyStr(sName.c_str());
                EdgeList[e].pname = CopyStr(cName.c_str());
                EdgeList[e].type = CS;
                CsN ++;

                e++;
                EdgeList[e].name = CopyStr(cName.c_str());
                EdgeList[e].pname = CopyStr(sName.c_str());

                EdgeList[e].node1 = (*it)->GetNodeIndex(2);
                EdgeList[e].node2 = (*it)->GetNodeIndex(3);
                EdgeList[e].lumped = 0;
                EdgeList[e].edgenum = e;
                EdgeList[e].type = CC;

                CcN++;
                EdgeList[e-1].pnum = e;
                EdgeList[e].pnum = e-1;
                break;
            case 'G':
                sName = (*it)->GetName();//"CS" + (*it)->GetName();
                cName = "VC_" + (*it)->GetName();  //(*it)->GetName();

		  	cout<<LINE_INFO<<"Element_G with CS name = "<<sName<<endl;
		  	cout<<LINE_INFO<<"Element_G with VC name = "<<cName<<endl;

                EdgeList[e].name = CopyStr(sName.c_str());
                EdgeList[e].pname = CopyStr(cName.c_str());
                EdgeList[e].type = CS;
                CsN ++;

                e++;
                EdgeList[e].name = CopyStr(cName.c_str());
                EdgeList[e].pname = CopyStr(sName.c_str());

                EdgeList[e].node1 = (*it)->GetNodeIndex(2);
                EdgeList[e].node2 = (*it)->GetNodeIndex(3);
                EdgeList[e].lumped = 0;
                EdgeList[e].edgenum = e;
                EdgeList[e].type = VC;

                VcN++;
                EdgeList[e-1].pnum = e;
                EdgeList[e].pnum = e-1;
                break;
            case 'H':
                sName = (*it)->GetName();//"VS" + (*it)->GetName();
                cName = (*it)->GetName();//"CC" + (*it)->GetName();

                EdgeList[e].name = CopyStr(sName.c_str());
                EdgeList[e].pname = CopyStr(cName.c_str());
                EdgeList[e].type = VS;
                VsN ++;

                e++;
                EdgeList[e].name = CopyStr(cName.c_str());
                EdgeList[e].pname = CopyStr(sName.c_str());

                EdgeList[e].node1 = (*it)->GetNodeIndex(2);
                EdgeList[e].node2 = (*it)->GetNodeIndex(3);
                EdgeList[e].lumped = 0;
                EdgeList[e].edgenum = e;
                EdgeList[e].type = CC;

                CcN++;
                EdgeList[e-1].pnum = e;
                EdgeList[e].pnum = e-1;
                break;
            default:
                cerr << "\nThe type of edge #" << e 
                    << " (" << EdgeList[e].name << ") "
                    << "incorrect.";
                delete [] EdgeList;
                EdgeList = 0;
                return 0;		// fail	
        }
	
		// Keep the max node num.
		if (max_n < EdgeList[e].node1)		max_n = EdgeList[e].node1;
		if (max_n < EdgeList[e].node2)		max_n = EdgeList[e].node2;

		++e;	// next edge
	}

	// Change by Xu Hui
	// Nodes are supposed to be numbered from 1 ... N //0 ... N-1.
	// So check whether (max_n > N). 
	if (max_n > N)
	{
	    cout << "\nNode number exceeded N-1. "
	        << "Check the node numbers in the graph file.";
	    delete [] EdgeList;
	    EdgeList = 0;
	    return 0;		// fail
	}

	cout << LINE_INFO << "End" << endl;
	
	return 1;
}
#endif

/** AddSymbol()
	adds a symbol to the head of the linked list  pointed by 
	(*SymbListHead) while checking the parallelism.
	Parallel symbols are lumped in this function by creating a 2D linked list.

	The leading parallel edge is given a name prefixed by "^".
	This function is called by GRDD::CreateSymbolList().

	The function also updates the counts: TotalNumSymbs; NumWorkingSymbs; NumWorkingEdges.

	(! Edge lumping is treated in this function.)
*/

void GPDD::AddSymbol(Edge *p_edge)
{
	// Always point at the tail of symblist. 
	// Use static to maintain the tail pointer.
	static SymbNode *pTail = NULL;  

	// First check whether this edge is in parallel to a symbol
	// already in the symblist.
	//cerr << "\nAdding edge - " << p_edge->name;
	
	SymbNode *pNode = SymbListTemp;
	SymbNode *pNewNode;

	// Only lumps Y/Z edges.
	if (p_edge->type == Y || p_edge->type == Z)
	{
		// Check the parallelism.
		while (pNode)
		{
	//cout << "\nFile: " << __FILE__ << ", Line: " << __LINE__;
			// Skip those non-Y/Z edges.
			if (pNode->edge->type > 1) 
			{
				pNode = pNode->next;  
				continue;
			}

			if ( (p_edge->node1 == pNode->edge->node1 && 
					p_edge->node2 == pNode->edge->node2) ||
					(p_edge->node1 == pNode->edge->node2 && 
					p_edge->node2 == pNode->edge->node1) )
			{
				//cerr << "\nDetected a parallel edge.";
				p_edge->lumped = 1;  	// This edge is lumped.

				// If lumping the first parallel edge, the
				// existing symbol is modified to be a composite symbol.
				if (!pNode->nextP)
				{
					//cerr << "\n !! The first lumping edge: " << p_edge->name;
					// Make a new node.
					pNewNode = new SymbNode;
					// Copy the 1st node to the 2nd.
					*pNewNode = *pNode;

					// Modify the symbol name of the composite symbol.
					char *new_name = new char[strlen(pNode->symbname)+2];
					sprintf(new_name, "^%s", pNode->symbname);
					pNode->symbname  = CopyStr(new_name);
					pNode->value = 0.0;		// unknown now
					// Should keep pNode->edge for parallelism checking.

					// Make connection.
					pNode->nextP = pNewNode;
					pNewNode->next  = NULL;
					pNewNode->nextP = NULL;
					delete[] new_name;
					new_name = 0;
				}

				// Go to the tail of the parallel symbol list.
				while (pNode->nextP)	pNode = pNode->nextP;
	//cout << "\nFile: " << __FILE__ << ", Line: " << __LINE__;
				// Add the new symbol node to the tail of parallel list.
				pNewNode = new SymbNode;
				pNode->nextP = pNewNode;

				pNewNode->symbname  = p_edge->name;
				pNewNode->value = p_edge->value;
				// cout << "value is : " << p_edge->value << endl;
				pNewNode->edge  = p_edge;
				pNewNode->next  = NULL;
				pNewNode->nextP = NULL;

				//cerr << "\nxxx A lumped edge added.";
				return;
			} // if 

			pNode = pNode->next;
		} // while
	}

	// Add the new symbol node to the tail of parallel list.
	++NumSymbols;		// Increase Nsymb
	
	// This edge is not lumped now.
	p_edge->lumped = 0;
	p_edge->index  = NumSymbols;		// Assign the index to this edge.

	pNewNode = new SymbNode;
	// Link to the symbol list.
	if (!SymbListTemp)	// The first symbol node.
		SymbListTemp = pTail = pNewNode;
	else			// The second on.
	{
		pTail->next = pNewNode,  pTail = pTail->next;  // Update the tail.
		//pNewNode->next = SymbListTemp;
		//SymbListTemp = pNewNode;
	}

	pNewNode->symbname  = p_edge->name;
	pNewNode->value = p_edge->value;
	pNewNode->edge  = p_edge;
	// cout << "value is : " << p_edge->value << endl;
	pNewNode->next  = NULL;
	pNewNode->nextP = NULL;
//	pNewNode->edgenum = edgeid;

	if(p_edge->type > 1)
	{
		pNewNode->isSrc = 1;
		pNewNode->pairedge = &EdgeList[pNewNode->edge->pnum];
	}
	else
	{
		pNewNode->isSrc = 0;
		pNewNode->pairedge = NULL;
	}
	//cerr << "\n--- New symbol node added.";
	return;
}


/**
  Look up the lumped symbol sublist of <tag> by matching the <name> and
  change the matched symbol value.
*/
void GPDD::ChangeSymbolValue(int tag, const char *name, double value)
{
//	cout << "value is : " << value << "tag is : " << tag << "name is : " << name <<endl;
	SymbNode *pNode;

	if (strcmp(LumpedSymbList[tag]->symbname, name) == 0)
	{
		LumpedSymbList[tag]->value = value;
	//	cout << "no lump.. " << endl;
	}
	else
	{
		pNode = LumpedSymbList[tag]->nextP;
		for (; pNode; pNode = pNode->nextP)
		{
//			cout << "initial value is : " << pNode->value << endl;
			if (strcmp(pNode->symbname, name) == 0)
			{
				pNode->value = value;
				break;
			}
		}
	}
		
	pNode = LumpedSymbList[tag]->nextP;
}

/**
	When a graph has been reduced into only one node, there could be some
	edges attached to that node as loops. Then we should check whether there
	exists any VS or CC edge. If yes, then the spanning tree (pair) so far created
	has CC or VS edge missing, hence cannot be valid.

	Rule: All CC and VS edges must appear in all spanning tree (pairs)
	        either in paired edges or common edges.
*/
int GPDD::CC_VS_Left(EdgeInfo *edgelist, int numedges)
{
	int 	num;
	
	for (int i = 0; i < numedges; i ++)
	{	
		num = edgelist[i].edgenum;
		
		if (EdgeList[num].type == VS || EdgeList[num].type == CC) {
			return ( YES);	// CC or VS in the tree ==> not valid
		}
	}

	return (NOPE);	// no CC and VS in tree ==> valid
}


/**	CollectNodes()
	collects all the GRDD nodes and put them together in the list ReduceList[] 
	for GRDD reduction.

	ReduceList[] holds the pointers of the GRDD nodes in a 2D linked list structure 
	such that those GRDD nodes of the same index "k" are linked as a sublist
	at ReduceList[k+2].
*/

void GPDD::CollectNodes()
{
	cout << LINE_INFO << endl;
	
///	UnMarkAllNodes(pGPDDroot, 1);
	
	ReduceList = new GPDDnode*[NumSymbols + 2];
	
	for (int i = 0; i < NumSymbols + 2; i ++)
		ReduceList[i] = NULL;
	
	ReduceList[0] = VertexZero;  // Put the two terminal nodes in the first two slots.
	ReduceList[1] = VertexOne;

	TraversalMark = pGPDDroot->GetMark();  // Get the root mark
	FlipTraversalMark();

	VertexOne->SetMark(TraversalMark);
	VertexZero->SetMark(TraversalMark);

	CollectNodes_R(pGPDDroot);

	pGPDDroot->SetMark(TraversalMark);  // The root has not been marked
}


/**
  An illustration of the ReduceList[]:
  ReduceList[0] ---[1] ----[2] ----[3] ---- .... [index] ...
                                       |          |                  |
                                       A1       B1                 F1
                                       A2       B2                 F2
                                       A3       B3                 F3
  where Ai's (resp. Bi's, Fi's) have the same symbol index, respectively.                                      
*/
void GPDD::CollectNodes_R(GPDDnode *p)
{
	if (!p || p->GetMark() == TraversalMark)
		return;

	p->SetMark(TraversalMark);
	p->SetNextNode(NULL);

	GPDDnode * q;

	// The next three lines add a GPDDnode to the sublist lead.
	q = ReduceList[p->GetSymbolIndex() + 2];
	p->SetNextNode(q);	
	ReduceList[p->GetSymbolIndex() + 2] = p;
	
	CollectNodes_R(p->GetVertex_IN());
	CollectNodes_R(p->GetVertex_OUT());
}

void GPDD::CountGPDDTerms()
{
	if(!pGPDDroot)
		NumTerms = 0;
	
//	UnMarkAllNodes(pGPDDroot, 1);

	VertexZero->SetNumTerms(0.0);
	VertexOne->SetNumTerms(1.0);

	TraversalMark = pGPDDroot->GetMark();
	FlipTraversalMark();
	
	VertexZero->SetMark(TraversalMark);
	VertexOne->SetMark(TraversalMark);

	NumTerms = CountGPDDTerms_R(pGPDDroot);

	pGPDDroot->SetMark(TraversalMark);
}

	
//int GPDD::CalculateTermsR(GPDDnode *p)
double GPDD::CountGPDDTerms_R(GPDDnode *p)
{
	if (!p)
		return 0;
	
	if (p->GetMark() == TraversalMark)
		return p->GetNumTerms();
	
	double terml, termr;
	terml = CountGPDDTerms_R(p->GetVertex_IN());
	termr = CountGPDDTerms_R(p->GetVertex_OUT());
	
	p->SetMark(TraversalMark);
	p->SetNumTerms(terml + termr);
	
	return (terml + termr);
}


/**
	A utility function for creating s-coefficients during s-expansion
	by successively updating the s-coefficient vector.
*/
//void GPDD::scoeffRNew(ExpGPDDnode * p)	// Old name by W. Chen
void GPDD::CreateCoeffList(ExpGPDDnode * p)
{
	vector<sCoeffNode*> *pCoeffShort	= p->Vertex_IN->sCoeffList;
	vector<sCoeffNode*> *pCoeffOpen = p->Vertex_OUT->sCoeffList;
	
	vector<int> *pPowerShort = p->Vertex_IN->sPowerList;
	vector<int> *pPowerOpen = p->Vertex_OUT->sPowerList;

	p->sCoeffList	= new vector<sCoeffNode*>;
	p->sPowerList = new vector<int>;

	int i, j, tag;

	int shortSize, openSize;

	shortSize	 = pPowerShort->size();
	openSize = pPowerOpen->size();

	sCoeffNode *pNewsCoeffNode, *pSharesCoeffNode;
	// sCoeffNode *pZNode = NULL;
	// sCoeffNode *pCNode = NULL;
	// sCoeffNode *pLNode = NULL;
	// ExpGPDDnode *pExpGPDDnode;

	char eType = LumpedSymbList[p->SymbolIndex]->symbname[0];
	
	if (eType == '^')
		eType = LumpedSymbList[p->SymbolIndex]->symbname[1];

	switch(eType)
	{
	case 'L':
		for (i = 0, j =0; i < shortSize || j < openSize;)
		{
			if (i < shortSize && j < openSize)
			{
				if ((*pPowerShort)[i] < ((*pPowerOpen)[j] + 1))
				{
					p->sPowerList->push_back(((*pPowerShort)[i] - 1));

					pNewsCoeffNode = new sCoeffNode;
					pNewsCoeffNode->pExpNode	= p;
					
					pNewsCoeffNode->Vertex_IN	= (*pCoeffShort)[i];
					pNewsCoeffNode->Vertex_OUT	= sCoeffZero;
					
					pSharesCoeffNode = Share_sCoeffNode(pNewsCoeffNode, tag);
					p->sCoeffList->push_back(pSharesCoeffNode);
					i ++;
				}
				else if ((*pPowerShort)[i] > ((*pPowerOpen)[j] + 1))
				{
					p->sPowerList->push_back((*pPowerOpen)[j]);

					pNewsCoeffNode = new sCoeffNode;
					pNewsCoeffNode->pExpNode	= p;
					
					pNewsCoeffNode->Vertex_IN	= sCoeffZero;
					pNewsCoeffNode->Vertex_OUT	= (*pCoeffOpen)[j];
					
					pSharesCoeffNode = Share_sCoeffNode(pNewsCoeffNode, tag);
					p->sCoeffList->push_back(pSharesCoeffNode);
					j ++;
				}
				else
				{
					p->sPowerList->push_back((*pPowerOpen)[j] );
					
					pNewsCoeffNode = new sCoeffNode;
					pNewsCoeffNode->pExpNode	= p;
					
					pNewsCoeffNode->Vertex_IN	= (*pCoeffShort)[i];
					pNewsCoeffNode->Vertex_OUT	= (*pCoeffOpen)[j];
					
					pSharesCoeffNode = Share_sCoeffNode(pNewsCoeffNode, tag);
					p->sCoeffList->push_back(pSharesCoeffNode);
					i ++;
					j ++;
				}
			}
			else if (i < shortSize)
			{
				p->sPowerList->push_back(((*pPowerShort)[i] - 1));
					
				pNewsCoeffNode = new sCoeffNode;
				pNewsCoeffNode->pExpNode	= p;
				pNewsCoeffNode->Vertex_OUT	= sCoeffZero;
				pNewsCoeffNode->Vertex_IN	= (*pCoeffShort)[i];
				pSharesCoeffNode = Share_sCoeffNode(pNewsCoeffNode, tag);
				p->sCoeffList->push_back(pSharesCoeffNode);
				i ++;
			}
			else
			{
				p->sPowerList->push_back((*pPowerOpen)[j]);
						
				pNewsCoeffNode = new sCoeffNode;
				pNewsCoeffNode->pExpNode	= p;
				pNewsCoeffNode->Vertex_OUT	= (*pCoeffOpen)[j];
				pNewsCoeffNode->Vertex_IN	= sCoeffZero;
				pSharesCoeffNode = Share_sCoeffNode(pNewsCoeffNode, tag);
				p->sCoeffList->push_back(pSharesCoeffNode);
				j ++;
			}
		}
		break;
		
	case 'C':
		for (i = 0, j =0; i < shortSize || j < openSize;)
		{
			if (i < shortSize && j < openSize)
			{
				if ((*pPowerShort)[i] < ((*pPowerOpen)[j] - 1))
				{
					p->sPowerList->push_back(((*pPowerShort)[i] + 1));

					pNewsCoeffNode = new sCoeffNode;
					pNewsCoeffNode->pExpNode	= p;
					
					pNewsCoeffNode->Vertex_IN	= (*pCoeffShort)[i];
					pNewsCoeffNode->Vertex_OUT	= sCoeffZero;
					
					pSharesCoeffNode = Share_sCoeffNode(pNewsCoeffNode, tag);
					p->sCoeffList->push_back(pSharesCoeffNode);
					i ++;
				}
				else if ((*pPowerShort)[i] > ((*pPowerOpen)[j] - 1))
				{
					p->sPowerList->push_back((*pPowerOpen)[j]);

					pNewsCoeffNode = new sCoeffNode;
					pNewsCoeffNode->pExpNode = p;
					
					pNewsCoeffNode->Vertex_IN	= sCoeffZero;
					pNewsCoeffNode->Vertex_OUT	= (*pCoeffOpen)[j];
					
					pSharesCoeffNode = Share_sCoeffNode(pNewsCoeffNode, tag);
					p->sCoeffList->push_back(pSharesCoeffNode);
					j ++;
				}
				else
				{
					p->sPowerList->push_back((*pPowerOpen)[j]);

					pNewsCoeffNode = new sCoeffNode;
					pNewsCoeffNode->pExpNode	= p;
					
					pNewsCoeffNode->Vertex_IN	= (*pCoeffShort)[i];
					pNewsCoeffNode->Vertex_OUT	= (*pCoeffOpen)[j];
					
					pSharesCoeffNode = Share_sCoeffNode(pNewsCoeffNode, tag);
					p->sCoeffList->push_back(pSharesCoeffNode);
					i ++;
					j ++;
				}
			}
			else if (i < shortSize)	// j = openSize
			{
				p->sPowerList->push_back(((*pPowerShort)[i] + 1));
				
				pNewsCoeffNode = new sCoeffNode;
				pNewsCoeffNode->pExpNode	= p;
				
				pNewsCoeffNode->Vertex_IN	= (*pCoeffShort)[i];
				pNewsCoeffNode->Vertex_OUT	= sCoeffZero;
				
				pSharesCoeffNode = Share_sCoeffNode(pNewsCoeffNode, tag);
				p->sCoeffList->push_back(pSharesCoeffNode);
				i ++;
			}
			else
			{
				p->sPowerList->push_back((*pPowerOpen)[j]);
					
				pNewsCoeffNode = new sCoeffNode;
				pNewsCoeffNode->pExpNode	= p;
				
				pNewsCoeffNode->Vertex_IN	= sCoeffZero;
				pNewsCoeffNode->Vertex_OUT	= (*pCoeffOpen)[j];
				
				pSharesCoeffNode = Share_sCoeffNode(pNewsCoeffNode, tag);
				p->sCoeffList->push_back(pSharesCoeffNode);
				j ++;
			}

		}
		break;
		
	default:  // R or G or other gains (independ on s)
		for(i = 0, j =0; i < shortSize || j < openSize;)
		{
			if (i < shortSize && j < openSize)
			{
				if ((*pPowerShort)[i] < (*pPowerOpen)[j])
				{
					p->sPowerList->push_back((*pPowerShort)[i]);

					pNewsCoeffNode = new sCoeffNode;
					pNewsCoeffNode->pExpNode	= p;
					pNewsCoeffNode->Vertex_OUT	= sCoeffZero;
					pNewsCoeffNode->Vertex_IN	= (*pCoeffShort)[i];
					pSharesCoeffNode = Share_sCoeffNode(pNewsCoeffNode, tag);
					p->sCoeffList->push_back(pSharesCoeffNode);
					i ++;
				}
				else if ((*pPowerShort)[i] > (*pPowerOpen)[j])
				{
					p->sPowerList->push_back((*pPowerOpen)[j]);

					pNewsCoeffNode = new sCoeffNode;
					pNewsCoeffNode->pExpNode	= p;
					pNewsCoeffNode->Vertex_OUT	= (*pCoeffOpen)[j];
					pNewsCoeffNode->Vertex_IN	= sCoeffZero;
					pSharesCoeffNode = Share_sCoeffNode(pNewsCoeffNode, tag);
					p->sCoeffList->push_back(pSharesCoeffNode);
					j ++;
				}
				else
				{
					p->sPowerList->push_back((*pPowerShort)[i]);

					pNewsCoeffNode = new sCoeffNode;
					pNewsCoeffNode->pExpNode	= p;
					pNewsCoeffNode->Vertex_OUT	= (*pCoeffOpen)[j];
					pNewsCoeffNode->Vertex_IN	= (*pCoeffShort)[i];
					pSharesCoeffNode = Share_sCoeffNode(pNewsCoeffNode, tag);
					p->sCoeffList->push_back(pSharesCoeffNode);
					i ++;
					j ++;
				}
			}
			else if (i < shortSize)
			{
				p->sPowerList->push_back((*pPowerShort)[i]);

				pNewsCoeffNode = new sCoeffNode;
				pNewsCoeffNode->pExpNode	= p;
				pNewsCoeffNode->Vertex_OUT	= sCoeffZero;
				pNewsCoeffNode->Vertex_IN	= (*pCoeffShort)[i];
				pSharesCoeffNode = Share_sCoeffNode(pNewsCoeffNode, tag);
				p->sCoeffList->push_back(pSharesCoeffNode);
				i ++;
			}
			else
			{
				p->sPowerList->push_back((*pPowerOpen)[j]);
				
				pNewsCoeffNode = new sCoeffNode;
				pNewsCoeffNode->pExpNode	= p;
				pNewsCoeffNode->Vertex_OUT	= (*pCoeffOpen)[j];
				pNewsCoeffNode->Vertex_IN	= sCoeffZero;
				pSharesCoeffNode = Share_sCoeffNode(pNewsCoeffNode, tag);
				p->sCoeffList->push_back(pSharesCoeffNode);
				j ++;
			}
		}
		break;
	}
}

int GPDD::CreateEdgeList(Edge_g * edge_list)
{
	cout << LINE_INFO << endl;
	
	if ( !edge_list )
	{
		cerr << LINE_INFO << "\nNo graph edge list available.";
		return   (STOP);
	}
	
	if (NumEdges > 0)
		EdgeList = new Edge[NumEdges];
	else
	{
		cerr << LINE_INFO << "\nError in CreateEdgeList(): No edge in graph.";
		return  (STOP);
	}
	
	Edge_g 	*edge;
	int		e = 0;
	int		max_n = 0;		// For checking the max node num.
	
	for (edge = edge_list; edge; edge = edge->next)
	{
		EdgeList[e].name = CopyStr( ToUpper(edge->name) ); 
		EdgeList[e].pname = NULL;
		EdgeList[e].node1 = edge->node1->num;
  		EdgeList[e].node2 = edge->node2->num;
		EdgeList[e].lumped = 0;	// Not lumped
		EdgeList[e].edgenum = e;
  		
		// Check syntax.
		if (EdgeList[e].node1 == EdgeList[e].node2)
		{
			cout << LINE_INFO << "\nError: Edge #" << e 
				 << " (" << EdgeList[e].name << ") "
				 << "is a self-loop.";
			delete [] EdgeList;
			EdgeList = 0;
			
			return (STOP);		// fail
		}

		// Make sure that the name of a pairing edge is set.
		if ( (edge->type >= CC) && (edge->pname == NULL) )
		{
			cerr << LINE_INFO << "\nThe name for one pairing edge is missing.";
		}
		EdgeList[e].type = edge->type;
		
		switch (edge->type)
		{
		case Z:
		case Y:
			EdgeList[e].value = edge->value;
			// cout << "edge value is : " << edge->value << endl;
			break;
		case VS:
			EdgeList[e].pname = CopyStr( ToUpper(edge->pname) );
			EdgeList[e].value = edge->value;
			VsN ++;
			break;
		case CS:
			EdgeList[e].pname = CopyStr( ToUpper(edge->pname) );
			EdgeList[e].value = edge->value;
			CsN ++;
			break;
		case VC:	
			EdgeList[e].pname = CopyStr( ToUpper(edge->pname) );
			VcN ++;
			break;
		case CC:
			EdgeList[e].pname = CopyStr( ToUpper(edge->pname) );
			CcN ++;
			break;
		default:
			cerr  << LINE_INFO << "\nThe type of edge #" << e 
				 << " (" << EdgeList[e].name << ") "
				 << "incorrect.";
			delete [] EdgeList;
			EdgeList = 0;
			
			return (STOP);		// fail	
		}		
		
		// Keep the max node num.
		max_n = (max_n < EdgeList[e].node1)	? EdgeList[e].node1 : max_n;
		max_n = (max_n < EdgeList[e].node2)	? EdgeList[e].node2 : max_n;
//		if (max_n < EdgeList[e].node1)		max_n = EdgeList[e].node1;
//		if (max_n < EdgeList[e].node2)		max_n = EdgeList[e].node2;
		
		++e;	// next edge
	}

	// Nodes are supposed to be numbered from 0 ... N-1.
	// So check whether (max_n > N-1). 
	if (max_n > NumNodes - 1)
	{
		cout << LINE_INFO << "\nNode number exceeded N-1. "
			 << "Check the node numbers in the graph file.";
		delete [] EdgeList;
		EdgeList = 0;
		
		return (STOP);		// fail
	}

	cout << LINE_INFO << "Set edge pairing ... " << endl;
	SetPairing();

	cout << LINE_INFO << "Print edge list ... " << endl;
	PrintFullEdgeList();
	
	return (OK);		// success
}


/**
	The function creates a child node for the current GPDD node.
	Naming convention: The successive node relation is described by
		"father" and "son". Hence, graph nodes, edges, etc. associated
		with the father Vertex are prefixed with "father", while those
		associated with a child vertex are prefixed with "son".
		
		Therefore, you'll see names like:  fatherXXXX and sonXXXX.

	The rules applied to the graph pair for "SYMBOL_IN" are:
	                          (L-graph)                  |  (R-graph)
		For VCVS (E):  Short VS                  |   Short VC (Open VS)
		For CCCS (F):  Short CS (Open CC)    |  Short CC
		For VCCS (G):  Short CS                   |  Short VC
		For CCVS (H):  Short VS (Open CC)    |  Short CC (Open VS)

	The rules applied to the graph pair for "SYMBOL_OUT" are:
	                            (L-graph)               |   (R-graph)
		For VCVS (E):   Short VS                |   Short VS (Open VC)
		For CCCS (F):	  Short CC (Open CS)  |   Short CC
		For VCCS (G):	  Open CS                 |   Open VC
		For CCVS (H):	  Short CC; Short VS   |  Short CC; Short VS
*/

//GPDDnode * GPDD::CreateGPDDnode(int SymbIndex, int optype, GPDDnode *fatherNode)
GPDDnode * GPDD::CreateGPDDnode(GPDDnode *fatherVertex, int optype)
{
#ifdef  TRACE	
	cout << LINE_INFO << (optype == SYMBOL_IN ? "(SYMBOL_IN)" : "(SYMBOL_OUT)") << endl;
#endif

	GPDDnode 	*pNewNode;
	Subgraph 	*subgraph_L, *subgraph_R, *hashgraph_L, *hashgraph_R;
	EdgeInfo 		*edgelist_L = NULL, *edgelist_R = NULL;

	int	tag = 0;
	int	signL = 1, signR = 1, sign;
	
	int	fatherEdges_L, fatherEdges_R;
	int	fatherNodes_L, fatherNodes_R;
	int	sonEdges_L, sonEdges_R;		// num edges
	int	sonNodes_L, sonNodes_R;	// num nodes

	int 	fatherIndex = fatherVertex->GetSymbolIndex();
	int 	sonIndex = fatherIndex + 1;		// Index increases from `father' to `son'

	SymbNode * fatherSymbol = LumpedSymbList[fatherIndex];
	char		*fatherName = LumpedSymbList[fatherIndex]->symbname;

//#ifdef TRACE
	string arrowname = (optype == SYMBOL_IN) ? "(left)" : "(right)";
//#endif

	if (fatherIndex >= NumSymbols) {
		cout << LINE_INFO << "\n\tSymbol index (#" << fatherIndex
			<< ") exceeds the number of symbols." << endl;
		return NULL;
	}

	// The potential son Node name. In case a terminal node is reached, this node is not created.
	char		*sonName = NULL;

	if (sonIndex < NumSymbols) {
		sonName = LumpedSymbList[sonIndex]->symbname;
	}

	fatherNodes_L = fatherVertex->GetNumNodes_L();
	fatherNodes_R = fatherVertex->GetNumNodes_R();
	
	fatherEdges_L = fatherVertex->GetNumEdges_L();
	fatherEdges_R = fatherVertex->GetNumEdges_R();

	/* 	Fix the remaining numbers of edges <Enum_L> and <Enum_R> for the L-graph and R-graph.
	
		<fatherSymbol> is associated with the source edge.
		If the edge pair involves CC, then the L-graph reduces 2 edges regardless of the operation.
		If the edge pair involves VS, then the R-graph reduces 2 edges regardless of the operation.
	*/
	if (fatherSymbol->isSrc && fatherSymbol->pairedge->type == CC)	//XS-CC
		sonEdges_L = fatherEdges_L - 2;		// The remaining num of edges for the L-graph
	else
		sonEdges_L = fatherEdges_L - 1;

	if (fatherSymbol->isSrc && fatherSymbol->edge->type == VS)		//VS-XC
		sonEdges_R = fatherEdges_R - 2;		// The remaining num of edges for the R-graph
	else
		sonEdges_R = fatherEdges_R - 1;

	// Fix the sign
	if (LumpedSymbList[fatherIndex]->edge->type == VS && optype == SYMBOL_IN)
		sign = -1;
	else
		sign = 1;


#ifdef TRACE
		cout << LINE_INFO << "NodeCreationCount = " << NodeCreationCount << endl; 
		cout << LINE_INFO << "\n\tFather \"" << fatherName << "\" (#" << fatherIndex
			<< ")--" << arrowname << "--> Son \""<< sonName << " (#" << sonIndex << ") " 
			<< "\n\tnumEdges_L  (father) = " << fatherEdges_L
			<< "\n\tnumEdges_R  (father) = " << fatherEdges_R
			<< "\n\tnumEdges_L  (son) = " << sonEdges_L
			<< "\n\tnumEdges_R  (son) = " << sonEdges_R
			<< endl;

		cout << LINE_INFO << "\n\tCheck the graph-pair at the father symbol \"" 
			<< fatherName << "\"" << endl;
		cout << LINE_INFO << "Edges of the L-graph: " << endl;
		//PrintPartialEdgeList(edgelist_L, sonEdges_L);
		PrintPartialEdgeList(fatherVertex->GetGraph_L(), fatherVertex->GetNumEdges_L());

		cout << LINE_INFO << "Edges of the R-graph: " << endl;
		//PrintPartialEdgeList(edgelist_R, sonEdges_R);
		PrintPartialEdgeList(fatherVertex->GetGraph_R(), fatherVertex->GetNumEdges_R());
#endif


	//Change by Xu Hui
	if (sonEdges_L == 0)	// No edge remains after the operation
	{
		if (optype == SYMBOL_OUT)	// Exclude the symbol
		{
#ifdef TRACE
			cout <<  LINE_INFO << "V[" << fatherName << "]--" << arrowname << "--> Zero!" << endl;
#endif
			return VertexZero;
		}
		else // (optype == SYMBOL_IN)
		{
			Edge 	*pCurEdge = LumpedSymbList[fatherIndex]->edge;
			int 		isSrc = LumpedSymbList[fatherIndex]->isSrc;
			int 		n1, n2, Sn1, Sn2, Cn1, Cn2;
			tag  = 1;

			if (!isSrc)
			{
				n1 = fatherVertex->GetGraph_L()[0].node1;
				n2 = fatherVertex->GetGraph_L()[0].node2;
				//cout<<LINE_INFO<<" n1 = "<<n1<<", n2 = "<<n2<<endl;
				if (n2 > n1)
					signL = -1;
			}
			sonNodes_L = 1;
		}
	}
	else
	{
		///cout << LINE_INFO << "(before) signL = " << signL << endl;
		
		edgelist_L = CreateLeftGraph(optype, fatherVertex, sonEdges_L, 
					fatherNodes_L, sonNodes_L, signL);
		// <sonNodes_L> and <signL> are modified upon return.
		
		///cout << LINE_INFO << "(after) signL = " << signL << endl;

		if (edgelist_L == NULL)
		{

#ifdef TRACE
			cout << LINE_INFO << "\n\tCreateLeftGraph() returns NULL; terminating at [VertexZero]"
				<< "\n\tV[" << fatherName << "]--" << arrowname << "--> Zero!" << endl;
#endif
			return VertexZero;
		}
		
		if (sonNodes_L > CountValidEdges(edgelist_L, sonEdges_L) + 1 || 
				!DetectTree(edgelist_L, sonEdges_L, sonNodes_L))
		{

			if (edgelist_L) {
				delete[] edgelist_L;
				edgelist_L = 0;
			}
			//		cout << LINE_INFO << "Code is here: return <Zero>" << endl << flush;
#ifdef TRACE
			cout << LINE_INFO << "\n\tLeft tree is not valid; terminating VertexZero!"
				<< "\n\tV[" << fatherName<< "]--" << arrowname << "--> Zero!" << endl;
#endif
			return VertexZero;
		}	

		if (sonNodes_L == 1)
		{
			if (CC_VS_Left(edgelist_L, sonEdges_L)== NOPE ) {
				tag = 1;	// Will take action after checking the R-graph
			}
			else
			{
#ifdef TRACE
				cout << LINE_INFO << "Print LEFT graph of Symbol[" << fatherName << "]" << endl;
				PrintPartialEdgeList(edgelist_L, sonEdges_L);

				cout << LINE_INFO << "\n\tLeft tree is not valid, terminating at VertexZero!"
					"\"" << fatherName<< "\"--" << arrowname << "--> Zero!" << endl;
#endif

				if (edgelist_L) {
					delete[] edgelist_L;
					edgelist_L = 0;
				}

				return VertexZero;
			}
		}
	}


	//Change by Xu Hui
	if (sonEdges_R == 0)	//only exist one edge or paired edges
	{
		if (optype == SYMBOL_OUT)	//Out the symbol
		{
#ifdef TRACE
			cout << LINE_INFO << "V[" << fatherName<< "]--" << arrowname << "--> Zero!" << endl;
#endif
			return VertexZero;
		}
		else
		{
			int isSrc = LumpedSymbList[fatherIndex]->isSrc;
			int n1, n2, Sn1, Sn2, Cn1, Cn2;

			if (!isSrc)
			{
				n1 = fatherVertex->GetGraph_R()[0].node1;
				n2 = fatherVertex->GetGraph_R()[0].node2;

				//cout<<LINE_INFO<<" n1 = "<<n1<<", n2 = "<<n2<<endl;
				if (n2 > n1)
					signR = -1;
			}

			sign = sign*signL*signR;
			fatherVertex->SetSign_IN(sign);

			if (tag == 1) 
				return VertexOne;
		}
	}
	else
	{
		///cout << LINE_INFO << "(before) signR = " << signR << endl;
		
		edgelist_R = CreateRightGraph(optype, fatherVertex, sonEdges_R, 
					fatherNodes_R, sonNodes_R, signR);
		// <sonNodes_R> and <signR> are modified upon return.
		
		///cout << LINE_INFO << "(after) signR = " << signR << endl;

#ifdef TRACE
			cout << LINE_INFO << "\n\tafter CreateRightGraph() : sonNodes_R = " << sonNodes_R << endl;
#endif

		if (edgelist_R == NULL)
		{
#ifdef TRACE
			cout << LINE_INFO << "\n\tCreateRightGraph() returns NULL; terminating at [VertexZero]"
				<< "\n\tV[" << fatherName << "]--" << arrowname << "--> Zero!" << endl;
#endif

			if (edgelist_L) {
				delete[] edgelist_L;
				edgelist_L = 0;
			}

			return VertexZero;
		}

		if (sonNodes_R > CountValidEdges(edgelist_R, sonEdges_R) + 1 || 
				!DetectTree(edgelist_R, sonEdges_R, sonNodes_R))
		{
#ifdef TRACE
			cout << LINE_INFO << "Print LEFT graph of Symbol \"" << fatherName << "\"" << endl;
			PrintPartialEdgeList(edgelist_L, sonEdges_L);

			cout << LINE_INFO << "Print RIGHT graph of Symbol \"" << fatherName << "\"" << endl;
			PrintPartialEdgeList(edgelist_R, sonEdges_R);

			cout << LINE_INFO << "\n\tRight tree is not valid, terminating at VertexZero!"
				<< "\n\t\"" << fatherName << "\"--" << arrowname << "--> Zero!" << endl;
#endif

			if (edgelist_L) {
				delete[] edgelist_L;
				edgelist_L = 0;
			}
			if (edgelist_R) {
				delete[] edgelist_R;
				edgelist_R = 0;
			}

			return VertexZero;
		}


		if (sonNodes_R == 1)
		{
			if (CC_VS_Left(edgelist_R, sonEdges_R) == NOPE && tag)	//tag = 1 means sonNodes_L = 1
			{
				sign = signL * signR * sign;

				if (sign == -1 && optype == SYMBOL_IN)
				{
					fatherVertex->SetSign_IN(sign);
				}
				else if (sign == -1 && optype == SYMBOL_OUT)
				{
					fatherVertex->SetSign_OUT(sign);
				}

				if (edgelist_L) {
					delete[] edgelist_L;
					edgelist_L = 0;
				}

				if (edgelist_R) {
///					PrintPartialEdgeList(edgelist_R, sonEdges_R);

					delete[] edgelist_R; 	//XXXX possible BUG here !!! XXXX
					edgelist_R = 0;
				}

#ifdef TRACE
				cout << LINE_INFO << "V[" << LumpedSymbList[fatherIndex]->symbname 
				<< "]--" << arrowname << "--> One!" << endl << flush;
#endif

				return VertexOne;
			}
			else
			{
#ifdef TRACE
				cout << LINE_INFO << " ... Print out LEFT graph of Symbol[" << fatherName << "]" << endl;
				PrintPartialEdgeList(edgelist_L, sonEdges_L);

				cout << LINE_INFO << " ... Print out RIGHT graph of Symbol[" << fatherName << "]" << endl;
				PrintPartialEdgeList(edgelist_R, sonEdges_R);

				cout << LINE_INFO << "Right graph not valid, Return VertexZero!" << endl;
				cout << LINE_INFO << "V[" << fatherName << "]--" << arrowname << "--> Zero!" << endl;
#endif

				if (edgelist_L) {
					delete[] edgelist_L;
					edgelist_L = 0;
				}

				if (edgelist_R) {
					delete[] edgelist_R;
					edgelist_R = 0;
				}

				return VertexZero;
			}
		}
	}

	// Now create a new GPDDnode which is neither of the terminal nodes.

	sign = signL * signR * sign;

	//Add NumNodes to Subgraph
	subgraph_L = new Subgraph(edgelist_L, sonEdges_L, sonNodes_L, signL, NumNodes);
	subgraph_R = new Subgraph(edgelist_R, sonEdges_R, sonNodes_R, signR, NumNodes);

#ifdef TRACE
	cout << LINE_INFO << "Finished reducing the graph pair." << endl;
	cout << LINE_INFO << "Print LEFT graph of Symbol[" 
		<< LumpedSymbList[sonIndex]->symbname << "]" << endl;
	subgraph_L->PrintGraph(EdgeList);
	
	cout << LINE_INFO << "Print RIGHT graph of Symbol[" 
		<< LumpedSymbList[sonIndex]->symbname << "]" << endl;
	subgraph_R->PrintGraph(EdgeList);

	cout << LINE_INFO << "\"" << fatherName << "\"--" << arrowname << "--> \"" 
		<< sonName << "\"" << endl;
#endif

	hashgraph_L  = ShareSubgraph(subgraph_L);
	hashgraph_R = ShareSubgraph(subgraph_R);

	pNewNode = new GPDDnode(LumpedSymbList[sonIndex], +1, +1, hashgraph_L, hashgraph_R);
	
	if (sign == -1 && optype == SYMBOL_IN)
	{
		fatherVertex->SetSign_IN(sign);
	}
	else if (sign == -1 && optype == SYMBOL_OUT)
	{
		fatherVertex->SetSign_OUT(sign);
	}
	
	return pNewNode;
}

/**
 	Create the Left Graph for the Child node.
 	optype = 0 --> SYMBOL_IN;
 	optype = 1 --> SYMBOL_OUT;
 	 <sonNodes> returns the number of nodes in the son graph.

	This function only applies operations to the Left-graph edges according
	to the symbol type and the operation type.

	Implementation details: 
		The edges from the father graph are duplicated to the son graph,
		with the edge(s) involved in the operation treated specially
		(either shorted or removed). In the case of shorting or collapsing
		an edge, the terminals nodes (n1, n2) are merged to one (retaining n1) 
		by which all other remaining edges connected to (n2) are renumbered to (n1).

		In the case of paired edges, the modification to the graph is a little more delicate,
		but the basic principle remains the same. See the code for the details.

	This is one of the kernel functions for graph reduction: 
		CreateLeftGraph() and CreateRightGraph().

	The rules applied to the Left graphare:
	                          (L-graph)                  |  (L-graph)
	                          "SYMBOL_IN"              |  "SYMBOL_OUT"
		For VCVS (E):  Short VS                  |   Short VS
		For CCCS (F):  Short CS (Open CC)    |  Open CS (Short CC)
		For VCCS (G):  Short CS                   |  Open CS 
		For CCVS (H):  Short VS (Open CC)    |  Short VS (Short CC)

	Operation "SYMBOL_IN" on the L-Graph:
		Short the Controlled edge (XS) meanwhile open CC if it exists.
	Operation "SYMBOL_OUT" on the L-Graph:
		Do case by case.
*/

EdgeInfo* GPDD::CreateLeftGraph(int optype, 
		GPDDnode *fatherVertex, int sonEdges, int fatherNodes, int &sonNodes, int &signL)
{
//	cout << LINE_INFO << endl;
	
	if (sonEdges == 0) {
		return NULL;
	}
	
	EdgeInfo *sonEdgelist_L, *pCurEdgeInfo;
	Edge 	*pCurEdge, *pCurPairEdge;

	int	SymbIndex = fatherVertex->GetSymbolIndex();
	int 	fatherEdges_L = fatherVertex->GetNumEdges_L();

/******************************************************************************/
//   S... : Controlled sources    Sn1: node1 of Controlled Sources   Sn2: node2 of Controlled Sources
//   C... : Controlling sources   Cn1: node1 of Controlling Sources   Cn2: node2 of Controlling Sources
//						    Se: Controlled source edges     Ce: Controlling source edges
/******************************************************************************/

	int 	isSrc = 0;
	int 	Stype = 0, Ctype = 0;
	int 	i = 0, j = 0;
	int 	tmp = 0;
	int	e = 0;
	int	n1 = 0, n2 = 0;
	int 	Sn1 = 0, Sn2 = 0;			// Terminal nodes of the Source edge
	int	Cn1 = 0, Cn2 = 0;			// Terminal nodes of the Controlling edge
	int	Se = 0, Ce = 0;			// numbers of Src & Ctrl edges
	
	//cout<<"NumEdge = "<<sonEdges<<endl;
	sonEdgelist_L = new EdgeInfo[sonEdges];
	
	pCurEdge 	= LumpedSymbList[SymbIndex]->edge;
	isSrc 		= LumpedSymbList[SymbIndex]->isSrc;
	
	signL = 1;
	sonNodes = fatherNodes;

/**
	First process the paired edges .................................................................
**/
	if (!isSrc)	//	Y/Z		
	{
		// Lookup the terminal node numbers (n1, n2) of the working edge.
		e = pCurEdge->edgenum;
		for (i = 0; i < fatherEdges_L; i ++)
		{
			if (fatherVertex->GetGraph_L()[i].edgenum == e)
				break;
		}
		n1 = fatherVertex->GetGraph_L()[i].node1;
		n2 = fatherVertex->GetGraph_L()[i].node2;

		if (optype == SYMBOL_IN)		
		{	
#ifdef  TRACE			
			cout << LINE_INFO << "\n\tWorking on a Y/Z edge for [SYMBOL_IN]" << endl;
#endif			
			if (n1 == n2)			// The working edge is already a loop
			{
				delete [] sonEdgelist_L;
				sonEdgelist_L = 0;
				return NULL;
			}
			
			if (n1 > n2)		// Always let (n1) < (n2); will retain the smaller node number.
			{
				tmp = n2,   n2 = n1,   n1 = tmp;
				signL = -signL;
			}

			/* Duplicate edges in the father graph to the son graph and do certain 
			    modifications like shorting edges for SYMBOL_IN op or opening edges
			    for SYMBOL_OUT op.
			*/
			for (j = 0, i = 0; i < fatherEdges_L; i ++)
			{
				if (j >= sonEdges) {
					cout << LINE_INFO << "Edges exceeds limit!" << endl;
					exit(0);
				}
				
				pCurEdgeInfo = &fatherVertex->GetGraph_L()[i];
				if (pCurEdgeInfo->edgenum == e)
				{
					// Find the edges for operation (short)
					if (GetPermNum(fatherVertex->GetGraph_L(), fatherEdges_L, n2))
					{
						signL = -signL;
					}
					continue;	// This edge is thus removed (not show up in the son graph)
				}

				// Enforce the effect of merging node (n2) to node (n1) and retaining (n1).
				sonEdgelist_L[j].node1 = (pCurEdgeInfo->node1 == n2) ? n1 : pCurEdgeInfo->node1;
				sonEdgelist_L[j].node2 = (pCurEdgeInfo->node2 == n2) ? n1 : pCurEdgeInfo->node2;

				sonEdgelist_L[j].edgenum = pCurEdgeInfo->edgenum;
				sonEdgelist_L[j].isSrc	   = pCurEdgeInfo->isSrc;
				sonEdgelist_L[j].mark = 0;

				j ++;
			}	
			sonNodes = fatherNodes -1;
		}
		else		// (optype == SYMBOL_OUT)
		{
#ifdef  TRACE			
			cout << LINE_INFO << "\n\tWorking on a Y/Z edge for [SYMBOL_OUT]" << endl;
#endif			
			// Duplicate the edges to the son graph, skipping the opened one.
			for (j = 0, i = 0; i < fatherEdges_L; i ++)
			{
				if (j >= sonEdges) {
					cout << LINE_INFO << "Edges exceeds limit!" << endl;
					exit(0);
				}

				pCurEdgeInfo = &fatherVertex->GetGraph_L()[i];
				if (pCurEdgeInfo->edgenum == e)
					continue;	// Skipping the edge
				
				sonEdgelist_L[j].node1 = pCurEdgeInfo->node1;
				sonEdgelist_L[j].node2 = pCurEdgeInfo->node2;
				
				sonEdgelist_L[j].edgenum = pCurEdgeInfo->edgenum;
				sonEdgelist_L[j].isSrc = pCurEdgeInfo->isSrc;
				sonEdgelist_L[j].mark = 0;

				j ++;
			}	
		}

		return sonEdgelist_L;
	}

/**
	Then process the paired edges .................................................................
**/
	// Preprocessing: To exchange nodes, get the operating edges, etc.
	
	pCurPairEdge = LumpedSymbList[SymbIndex]->pairedge;
	Stype = pCurEdge->type;
	Ctype = pCurPairEdge->type;
	
	Se = pCurEdge->edgenum;
	Ce = pCurPairEdge->edgenum; 
	/** Only the CC edge from the R-graph is allowed in the L-graph. **/


#ifdef TRACE
		cout << LINE_INFO << "\n\tSe = " << Se << ", Ce = " << Ce << endl;
#endif


	/** Look up the Source edge index. **/
	for (i = 0; i < fatherEdges_L; i ++)
	{
		if (fatherVertex->GetGraph_L()[i].edgenum == Se)
			break;
	}

	if (i == fatherEdges_L) {
		cout << LINE_INFO << "\n\tFather edge (XS) was not found in L-graph." 
			<< "\n\tProgram is terminated." << endl;
		exit(0);
	}
	
	Sn1 = fatherVertex->GetGraph_L()[i].node1;  
	Sn2 = fatherVertex->GetGraph_L()[i].node2;

	/** An VS edge must be included in any admissible tree.
		The VS edge is "shorted" regardless of the "SYMBOL_IN" or "SYMBOL_OUT"
		operation. Hence, a VS edge is not allowed to becom a loop.
	**/
	if (Stype == VS && Sn1 == Sn2)
	{
		delete[] sonEdgelist_L;
		sonEdgelist_L = 0;
		return NULL;			// invalid tree
	}
	
	/** Look up the CC edge which must be included in the L-graph as well. **/
	if (Ctype == CC)
	{
		for (i = 0; i < fatherEdges_L; i ++)
		{
			if (fatherVertex->GetGraph_L()[i].edgenum == Ce)
				break;
		}

		if (i == fatherEdges_L) {
			cout << LINE_INFO << "\n\tFather edge (CC) was not found in L-graph."
				<< "\n\tProgram is terminated." << endl;
			exit(0);
		}

		Cn1 = fatherVertex->GetGraph_L()[i].node1;
		Cn2 = fatherVertex->GetGraph_L()[i].node2;
	} else {
		Cn1 = Cn2 = -1;	//  Cn1 and Cn2 are used only in the case that CC_S.
	}

	//*********** Work on the operations for the dependent edge pair ************8
	if (optype == SYMBOL_IN)				
	{
#ifdef  TRACE			
			cout << LINE_INFO << "\n\tWorking on a Controlled edge for [SYMBOL_IN]" << endl;
#endif			

		/**Whether VS is a loop was checked before. **/
		if (Stype == CS && Sn1 == Sn2)	  // The working source edge is a loop.
		{
			delete []sonEdgelist_L;
			sonEdgelist_L = 0;
			return NULL;			// Invalid tree
		}

		// The source edge is to be shorted
		if (Sn1 > Sn2)		// Let (Sn1) < Sn2)
		{
			tmp = Sn2,    Sn2 = Sn1,    Sn1 = tmp;
			signL = -signL;
		}

		for (j = 0, i = 0; i < fatherEdges_L; i ++)
		{
			if (j >= sonEdges) {
				cout << LINE_INFO << "Edges exceeds limit!" << endl;
				cout << "\tfatherEdges_L = " << fatherEdges_L 
					<< "; sonEdges = " << sonEdges << " (j = " << j << ")" << endl;
				exit(0);
			}

			pCurEdgeInfo = &fatherVertex->GetGraph_L()[i];


#ifdef TRACE
				cout << LINE_INFO << "\n\tSe = " << Se << ", Ce = " << Ce << endl;
				cout << "\tedgenum (info) = " << pCurEdgeInfo->edgenum << endl;
				cout << "\tedgenum = " << pCurEdge->edgenum << endl;
#endif

			
			if (pCurEdgeInfo->edgenum == Se || pCurEdgeInfo->edgenum == Ce)
			{
				if (pCurEdgeInfo->edgenum == Se && 
					(GetPermNum(fatherVertex->GetGraph_L(), fatherEdges_L, Sn2)))
				{
					signL = -signL;
				}
				continue;
				/** Skip the "shorted" source edge (CS/VS) and 
				      skip the "opened"  CC edge (if exists). **/
			}

			// Enforce the effect of merging node (Sn2) to node (Sn1) and retaining (Sn1).
			sonEdgelist_L[j].node1 = (pCurEdgeInfo->node1 == Sn2) ? Sn1 : pCurEdgeInfo->node1;
			sonEdgelist_L[j].node2 = (pCurEdgeInfo->node2 == Sn2) ? Sn1 : pCurEdgeInfo->node2;

			sonEdgelist_L[j].edgenum = pCurEdgeInfo->edgenum;
			sonEdgelist_L[j].isSrc = pCurEdgeInfo->isSrc;
			sonEdgelist_L[j].mark = 0;
			j ++;
		}

		/** The source edge (CS/VS) is collapsed, hence #nodes is reduced by 1.
		**/
		sonNodes = fatherNodes - 1;
		return sonEdgelist_L;
	}
	
	if (optype == SYMBOL_OUT)	// Implement the all four cases
	{
		/** CC is to be shorted. So it is already a loop, we should stop.
		**/
		if (Ctype == CC && Cn1 == Cn2)
		{
			delete[] sonEdgelist_L;
			sonEdgelist_L = 0;
			return NULL;			// Invalid tree.
		}

		if (Stype == VS && Ctype == VC)	// (E) : Short VS
		{
#ifdef  TRACE			
			cout << LINE_INFO << "\n\tWorking on a VC-VS edge for [SYMBOL_OUT]" << endl;
#endif			

			// The source edge is to be shorted
			if (Sn1 > Sn2)		// Let (Sn1) < Sn2)
			{
				tmp = Sn2,    Sn2 = Sn1,    Sn1 = tmp;
				signL = -signL;
			}

			for (j = 0, i = 0; i < fatherEdges_L; i ++)
			{
				if (j >= sonEdges) {
					cout << LINE_INFO << "Edges exceeds limit!" << endl;
					exit(0);
				}

				pCurEdgeInfo = &fatherVertex->GetGraph_L()[i];
				if (pCurEdgeInfo->edgenum == Se)
				{
					if (GetPermNum(fatherVertex->GetGraph_L(), fatherEdges_L, Sn2))
					{
						signL = -signL;
					}
					continue;	// Short the VS edge
				}

				// Enforce the effect of merging node (Sn2) to node (Sn1) and retaining (Sn1).
				sonEdgelist_L[j].node1 = (pCurEdgeInfo->node1 == Sn2) ? Sn1 : pCurEdgeInfo->node1;
				sonEdgelist_L[j].node2 = (pCurEdgeInfo->node2 == Sn2) ? Sn1 : pCurEdgeInfo->node2;

				sonEdgelist_L[j].edgenum = pCurEdgeInfo->edgenum;
				sonEdgelist_L[j].isSrc = pCurEdgeInfo->isSrc;
				sonEdgelist_L[j].mark = 0;
				j ++;

			}

			// The VS edge has been collapsed, hence #nodes reduces by one. 
			sonNodes = fatherNodes - 1;
			return sonEdgelist_L;
		}
		else if (Stype == CS && Ctype == CC)	// (F) : Open CS (Short CC)
		{
#ifdef  TRACE			
			cout << LINE_INFO << "\n\tWorking on a CC-CS edge for [SYMBOL_OUT]" << endl;
#endif			

			if (Cn1 > Cn2)   // force (Cn1) < (Cn2)
			{
				tmp = Cn2,    Cn2 = Cn1,    Cn1 = tmp;
				signL = -signL;
			}

			for (j = 0, i = 0; i < fatherEdges_L; i ++)
			{
				if (j >= sonEdges) {
					cout << LINE_INFO << "Edges exceeds limit!" << endl;
					exit(0);
				}

				pCurEdgeInfo = &fatherVertex->GetGraph_L()[i];
				if (pCurEdgeInfo->edgenum == Se || pCurEdgeInfo->edgenum == Ce)
				{
					if (pCurEdgeInfo->edgenum == Ce && 
						GetPermNum(fatherVertex->GetGraph_L(), fatherEdges_L, Cn2))
						signL = -signL;
					continue;	// Skip CS and CC
				}

				// Because CC is shorted,
				// enforce the effect of merging node (Cn2) to node (Cn1) and retaining (Cn1).
				sonEdgelist_L[j].node1 = (pCurEdgeInfo->node1 == Cn2) ? Cn1 : pCurEdgeInfo->node1;
				sonEdgelist_L[j].node2 = (pCurEdgeInfo->node2 == Cn2) ? Cn1 : pCurEdgeInfo->node2;

				sonEdgelist_L[j].edgenum = pCurEdgeInfo->edgenum;
				sonEdgelist_L[j].isSrc = pCurEdgeInfo->isSrc;
				sonEdgelist_L[j].mark = 0;
				j ++;

			}
			sonNodes = fatherNodes - 1;	// because CC is collapsed
			return sonEdgelist_L;
		}
		else if  (Stype == CS && Ctype == VC)	// (G) : Open CS
		{
#ifdef  TRACE			
			cout << LINE_INFO << "\n\tWorking on a VC-CS edge for [SYMBOL_OUT]" << endl;
#endif			

			for (j = 0, i = 0; i < fatherEdges_L; i ++)
			{
				if (j >= sonEdges) {
					cout << LINE_INFO << "Edges exceeds limit!" << endl;
					exit(0);
				}

				pCurEdgeInfo = &fatherVertex->GetGraph_L()[i];
				if (pCurEdgeInfo->edgenum == Se)
					continue;
				
				sonEdgelist_L[j].node1 = pCurEdgeInfo->node1;
				sonEdgelist_L[j].node2 = pCurEdgeInfo->node2;
				
				sonEdgelist_L[j].edgenum = pCurEdgeInfo->edgenum;
				sonEdgelist_L[j].isSrc = pCurEdgeInfo->isSrc;
				sonEdgelist_L[j].mark = 0;
				j ++;
			}
			return sonEdgelist_L;
		}
		else if  (Stype == VS && Ctype == CC)	// (H) : Short VS (Short CC)
		{
#ifdef  TRACE			
			cout << LINE_INFO << "\n\tWorking on a CC-VS edge for [SYMBOL_OUT]" << endl;
#endif			

			if (Cn1 == Sn1 && Cn2 == Sn2)   // The CC edge and VS edge together form a loop
			{
				delete[] sonEdgelist_L;
				sonEdgelist_L = 0;
				return NULL;
			}
			
			// The controlling edge is to be shorted
			if (Cn1 > Cn2)
			{
				tmp = Cn2,    Cn2 = Cn1,    Cn1 = tmp;
				signL = -signL;
			}

			// The source edge is to be shorted
			if (Sn1 > Sn2)		// Let (Sn1) < Sn2)
			{
				tmp = Sn2,    Sn2 = Sn1,    Sn1 = tmp;
				signL = -signL;
			}

			/* In case the two edges (CC & VS) are connected, consider the following cases:
			 //      Case 1:      (Sn1) --CS--(Sn2 | Cn2) --CC-- (Cn1)
			 //      Case 2:      (Sn1) --CS--(Sn2 | Cn1) --CC-- (Cn2)
			 //      Case 3:      (Sn2) --CS--(Sn1 | Cn2) --CC-- (Cn1)
			 //      Case 4:      (Sn2) --CS--(Sn1 | Cn1) --CC-- (Cn2)
			 
			    Note that (Sn1) < (Sn2) and (Cn1) < (Cn2) have been enforced. 
			    After arrangement, the exterior terminals are n1 and n2 with (n1) < (n2).
			    Then, renumber both (Sn1) & (Sn2) to (n1); renumber both (Cn1) & (Cn2) to (n2).
			**/
			if (Sn2 == Cn2)
			{
			      //   Case 1:      (Sn1) --CS--[Sn2 | Cn2] --CC-- (Cn1) ==> (Sn1) ----(Cn1)
				if (Sn1 < Cn1)
					n1 = Sn1,    n2 = Sn1;  // ==> (n1 == n2)
				else
					n1 = Cn1,    n2 = Cn1;	 // ==> (n1 == n2)
					//n2 = Cn1;		//Why not n2=Cn1?	(* corrected *)
			}
			else if (Sn2 == Cn1)
			       //  Case 2:      (Sn1) --CS--[Sn2 | Cn1] --CC-- (Cn2) ==> (Sn1) ----(Cn2)
				n1 = Sn1,    n2 = Sn1;	// ==> (n1 == n2)
			else if (Sn1 == Cn2)
			       //  Case 3:      (Sn2) --CS--[Sn1 | Cn2] --CC-- (Cn1) ==> (Sn2) ----(Cn1)
				n1 = Cn1,    n2 = Cn1;  // ==> (n1 == n2)
			else // For Case 4 and for the case the two edges are not concatenated.
			       //  Case 4:      (Sn2) --CS--[Sn1 | Cn1] --CC-- (Cn2) ==> (Sn2) ----(Cn2)
				n1 = Sn1,    n2 = Cn1;		// retain the smaller numbers

			for (j = 0, i = 0; i < fatherEdges_L; i ++)
			{
				if (j >= sonEdges) {
					cout << LINE_INFO << "Edges exceeds limit!" << endl;
					exit(0);
				}

				pCurEdgeInfo = &fatherVertex->GetGraph_L()[i];
				if (pCurEdgeInfo->edgenum == Se) 
				{
					if (GetPermNum(fatherVertex->GetGraph_L(), fatherEdges_L, Sn2))
						signL = -signL;
					continue;
				}
				else if (pCurEdgeInfo->edgenum == Ce)
				{
					if (GetPermNum(fatherVertex->GetGraph_L(), fatherEdges_L, Cn2))
						signL = -signL;
					continue;
				}


				// Enforce the effect of merging node (Sn2) to node (Sn1) and retaining (n1).
				sonEdgelist_L[j].node1 = (pCurEdgeInfo->node1 == Sn1 || pCurEdgeInfo->node1 == Sn2) ? 
							n1 : pCurEdgeInfo->node1;
				sonEdgelist_L[j].node2 = (pCurEdgeInfo->node2 == Sn1 || pCurEdgeInfo->node2 == Sn2) ? 
							n1 : pCurEdgeInfo->node2;

				sonEdgelist_L[j].edgenum = pCurEdgeInfo->edgenum;
				sonEdgelist_L[j].isSrc = pCurEdgeInfo->isSrc;
				sonEdgelist_L[j].mark = 0;
				j ++;

			}

			// Renumber both nodes (Cn1) and (Cn2) to (n2) by scanning the edges once again.
			for (i = 0; i < sonEdges; i ++)
			{
				pCurEdgeInfo = &sonEdgelist_L[i];
				
				if (pCurEdgeInfo->node1 == Cn1 || pCurEdgeInfo->node1 == Cn2)
					pCurEdgeInfo->node1 = n2;
				
				if (pCurEdgeInfo->node2 == Cn1 || pCurEdgeInfo->node2 == Cn2)
					pCurEdgeInfo->node2 =  n2;
			}
			sonNodes = fatherNodes - 2;  // Shorting two edges results in the reduction of 2 nodes
			return sonEdgelist_L;
		}
		else
		{
			delete[] sonEdgelist_L;
			sonEdgelist_L = 0;
			return NULL;
		}
	}
	
}

/**
 	Create the Left Graph for the Child node.
 	optype = 0 --> SYMBOL_IN;
 	optype = 1 --> SYMBOL_OUT;
 	 <sonNodes> returns the number of nodes in the son graph.

	This function only applies operations to the Right-graph edges according
	to the symbol type and the operation type.

	Implementation details: 
		The edges from the father graph are duplicated to the son graph,
		with the edge(s) involved in the operation treated specially
		(either shorted or removed). In the case of shorting or collapsing
		an edge, the terminals nodes (n1, n2) are merged to one (retaining n1) 
		by which all other remaining edges connected to (n2) are renumbered to (n1).

		In the case of paired edges, the modification to the graph is a little more delicate,
		but the basic principle remains the same. See the code for the details.

	This is one of the kernel functions for graph reduction: 
		CreateLeftGraph() and CreateRightGraph().

	The rules applied to the Right graph pair:
	                   (R-graph)                         |  (R-graph)
	                   "SYMBOL_IN"                     | "SYMBOL_OUT"
		For VCVS (E):  Short VC (Open VS)   |  Open VC (Short VS)
		For CCCS (F):  Short CC                  |  Short CC
		For VCCS (G):  Short VC                  |  Open VC
		For CCVS (H):  Short CC (Open VS)   |  Short CC (Short VS)

	Operation "SYMBOL_IN" on the R-Graph:
		Short the Controlling edge (XC) meanwhile open VS if it exists.
	Operation "SYMBOL_OUT" on the R-Graph:
		Do case by case.
*/
EdgeInfo* GPDD::CreateRightGraph(int optype, GPDDnode *fatherVertex, 
		int sonEdges, int fatherNodes, int &sonNodes, int& signR)
{
//	cout << LINE_INFO << endl;
	
	if (sonEdges == 0) {
		return NULL;
	}

	EdgeInfo *sonEdgelist_R, *pCurEdgeInfo;
	Edge 	*pCurEdge, *pCurPairEdge;

	int	SymbIndex = fatherVertex->GetSymbolIndex();
	int 	fatherEdges_R = fatherVertex->GetNumEdges_R();

/******************************************************************************/
//   S... : Controlled sources    Sn1: node1 of Controlled Sources   Sn2: node2 of Controlled Sources
//   C... : Controlling sources   Cn1: node1 of Controlling Sources   Cn2: node2 of Controlling Sources
//						    Se: Controlled source edges     Ce: Controlling source edges
/******************************************************************************/

	int 	isSrc = 0;
	int 	Stype = 0, Ctype = 0;
	int 	i = 0, j = 0;
	int 	tmp = 0;
	int	e = 0;
	int	n1 = 0, n2 = 0;
	int 	Sn1 = 0, Sn2 = 0;			// Terminal nodes of the Source edge
	int	Cn1 = 0, Cn2 = 0;			// Terminal nodes of the Controlling edge
	int	Se = 0, Ce = 0;			// numbers of Src & Ctrl edges
	
	sonEdgelist_R = new EdgeInfo[sonEdges];
	
	pCurEdge = LumpedSymbList[SymbIndex]->edge;
	isSrc 	= LumpedSymbList[SymbIndex]->isSrc;

	signR = 1;
	sonNodes = fatherNodes;

/**
	First process the paired edges .................................................................
**/
	if (!isSrc)
	{
		// Lookup the terminal node numbers (n1, n2) of the working edge.
		e = pCurEdge->edgenum;
		for (i = 0; i < fatherEdges_R; i ++)
		{
			if (fatherVertex->GetGraph_R()[i].edgenum == e)
				break;
		}
		n1 = fatherVertex->GetGraph_R()[i].node1;  //pCurEdge->node1;
		n2 = fatherVertex->GetGraph_R()[i].node2;

		if (optype == SYMBOL_IN)	
		{	
#ifdef  TRACE			
			cout << LINE_INFO << "\n\tWorking on a Y/Z edge for [SYMBOL_IN]" << endl;
#endif			

			if (n1 == n2)			// Already a loop
			{
				delete []sonEdgelist_R;
				sonEdgelist_R = 0;
				return NULL;
			}
			if (n1 > n2)
			{
				tmp = n2,   n2 = n1,   n1 = tmp;
				signR = -signR;
			}

			for (j = 0, i = 0; i < fatherEdges_R; i ++)
			{
				if (j >= sonEdges) {
					cout << LINE_INFO << "Edges exceeds limit!" << endl;
					exit(0);
				}

				pCurEdgeInfo = &fatherVertex->GetGraph_R()[i];
				if (pCurEdgeInfo->edgenum == e)
				{
					if (GetPermNum(fatherVertex->GetGraph_R(), fatherEdges_R, n2))
						signR = -signR;
					continue;
				}

				// Enforce the effect of merging node (n2) to node (n1) and retaining (n1).
				sonEdgelist_R[j].node1 = (pCurEdgeInfo->node1 == n2) ? n1 : pCurEdgeInfo->node1;
				sonEdgelist_R[j].node2 = (pCurEdgeInfo->node2 == n2) ? n1 : pCurEdgeInfo->node2;
				
				sonEdgelist_R[j].edgenum = pCurEdgeInfo->edgenum;
				sonEdgelist_R[j].isSrc = pCurEdgeInfo->isSrc;
				sonEdgelist_R[j].mark = 0;

				j ++;
			}	
			
			sonNodes = fatherNodes - 1;
		}
		else		// (optype = SYMBOL_OUT)
		{
#ifdef  TRACE			
			cout << LINE_INFO << "\n\tWorking on a Y/Z edge for [SYMBOL_OUT]" << endl;
#endif			

			// Duplicate the edges to the son graph, skipping the opened one.
			for (j = 0, i = 0; i < fatherEdges_R; i ++)
			{
				if (j >= sonEdges) {
					cout << LINE_INFO << "Edges exceeds limit!" << endl;
					exit(0);
				}

				pCurEdgeInfo = &fatherVertex->GetGraph_R()[i];
				if (pCurEdgeInfo->edgenum == e)
					continue;	// Skipping the edge
				
				sonEdgelist_R[j].node1 = pCurEdgeInfo->node1;
				sonEdgelist_R[j].node2 = pCurEdgeInfo->node2;
				
				sonEdgelist_R[j].edgenum = pCurEdgeInfo->edgenum;
				sonEdgelist_R[j].isSrc = pCurEdgeInfo->isSrc;
				sonEdgelist_R[j].mark = 0;

				j ++;
			}	
		}

		return sonEdgelist_R;
	}


/**
	Then process the paired edges .................................................................
**/
	// Preprocessing: To exchange nodes, get the operating edges, etc.

	pCurPairEdge = LumpedSymbList[SymbIndex]->pairedge;
	Stype = pCurEdge->type;
	Ctype = pCurPairEdge->type;
	
	Se = pCurEdge->edgenum;
	Ce = pCurPairEdge->edgenum; 
	
	for (i = 0; i < fatherEdges_R; i ++)
	{
		if (fatherVertex->GetGraph_R()[i].edgenum == Ce)
			break;
	}

	if (i == fatherEdges_R) {
		cout << LINE_INFO << "\n\tFather edge (XC) was not found in R-graph."
			<< "\n\tProgram is terminated." << endl;
		exit(0);
	}

	Cn1 = fatherVertex->GetGraph_R()[i].node1;
	Cn2 = fatherVertex->GetGraph_R()[i].node2;

	/** A CC edge must be included in any admissible tree.
		The CC edge is "shorted" regardless of the "SYMBOL_IN" or "SYMBOL_OUT"
		operation. Hence, a CC edge is not allowed to becom a loop.
	**/
	if (Ctype == CC && Cn1 == Cn2)		// This edge is a loop
	{
		delete []sonEdgelist_R;
		sonEdgelist_R = 0;
		return NULL;
	}
	
	/** Look up the VS edge which must be included in the R-graph as well. **/
	if (Stype == VS)
	{
		for (i = 0; i < fatherEdges_R; i ++)
		{
			if (fatherVertex->GetGraph_R()[i].edgenum == Se)
				break;
		}

		if (i == fatherEdges_R) {
			cout << LINE_INFO << "\n\tFather edge (VS) was not found in R-graph."
				<< "\n\tProgram is terminated." << endl;
			exit(0);
		}
		
		Sn1 = fatherVertex->GetGraph_R()[i].node1;  //pCurEdge->node1;
		Sn2 = fatherVertex->GetGraph_R()[i].node2;
	} else {
		Sn1 = Sn2 = -1;
	}

	//*********** Work on the operations for the dependent edge pair ************8
	if (optype == SYMBOL_IN)				// select this symbol
	{
#ifdef  TRACE			
			cout << LINE_INFO << "\n\tWorking on a Controlling edge for [SYMBOL_IN]" << endl;
#endif			

		if (Ctype == VC && Cn1 == Cn2)		// This is a loop
		{
			delete []sonEdgelist_R;
			sonEdgelist_R = 0;
			return NULL;
		}		

		// The controlling edge is to be shorted
		if (Cn1 > Cn2)
		{
			tmp = Cn2,    Cn2 = Cn1,    Cn1 = tmp;
			signR = -signR;
		}

		for (j = 0, i = 0; i < fatherEdges_R; i ++)
		{
			if (j >= sonEdges) {
				cout << LINE_INFO << "Edges exceeds limit!" << endl;
				exit(0);
			}

			pCurEdgeInfo = &fatherVertex->GetGraph_R()[i];
			if (pCurEdgeInfo->edgenum == Se || pCurEdgeInfo->edgenum == Ce)
			{
				if (GetPermNum(fatherVertex->GetGraph_R(), fatherEdges_R, Cn2) 
						&& pCurEdgeInfo->edgenum == Ce)
				{
					signR = -signR;
				}
				continue;
			}

			// Enforce the effect of merging node (n2) to node (n1) and retaining (n1).
			sonEdgelist_R[j].node1 = (pCurEdgeInfo->node1 == Cn2) ? Cn1 : pCurEdgeInfo->node1;
			sonEdgelist_R[j].node2 = (pCurEdgeInfo->node2 == Cn2) ? Cn1 : pCurEdgeInfo->node2;

			sonEdgelist_R[j].edgenum = pCurEdgeInfo->edgenum;
			sonEdgelist_R[j].isSrc = pCurEdgeInfo->isSrc;
			sonEdgelist_R[j].mark = 0;
			j ++;
		}
		sonNodes = fatherNodes - 1;
		return sonEdgelist_R;
	}
	
	if (optype == SYMBOL_OUT)
	{
		if (Stype == VS && Sn1 == Sn2)		// already a loop
		{
			delete[] sonEdgelist_R;
			sonEdgelist_R = 0;
			return NULL;
		}
		
		if (Stype == VS && Ctype == VC) // (E) : short VS, open VC
		{
#ifdef  TRACE			
			cout << LINE_INFO << "\n\tWorking on a VC-VS edge for [SYMBOL_OUT]" << endl;
#endif			

			if (Sn1 > Sn2)
			{
				tmp = Sn2,    Sn2 = Sn1,    Sn1 = tmp;
				signR = -signR;
			}

			for (j = 0, i = 0; i < fatherEdges_R; i ++)
			{
				if (j >= sonEdges) {
					cout << LINE_INFO << "Edges exceeds limit!" << endl;
					exit(0);
				}

				pCurEdgeInfo = &fatherVertex->GetGraph_R()[i];
				if (pCurEdgeInfo->edgenum == Se || pCurEdgeInfo->edgenum == Ce)
				{
					if (GetPermNum(fatherVertex->GetGraph_R(), fatherEdges_R, Sn2) && pCurEdgeInfo->edgenum == Se)
						signR = -signR;
					continue;
				}

				// Enforce the effect of merging node (n2) to node (n1) and retaining (n1).
				sonEdgelist_R[j].node1 = (pCurEdgeInfo->node1 == Sn2) ? Sn1 : pCurEdgeInfo->node1;
				sonEdgelist_R[j].node2 = (pCurEdgeInfo->node2 == Sn2) ? Sn1 : pCurEdgeInfo->node2;

				sonEdgelist_R[j].edgenum = pCurEdgeInfo->edgenum;
				sonEdgelist_R[j].isSrc = pCurEdgeInfo->isSrc;
				sonEdgelist_R[j].mark = 0;
				j ++;

			}
			sonNodes = fatherNodes - 1;
			return sonEdgelist_R;
		}
		else if (Stype == CS && Ctype == CC)  // (F) : open CS, short CC
		{
#ifdef  TRACE			
			cout << LINE_INFO << "\n\tWorking on a CC-CS edge for [SYMBOL_OUT]" << endl;
#endif			

			// The controlling edge is to be shorted
			if (Cn1 > Cn2)
			{
				tmp = Cn2,    Cn2 = Cn1,    Cn1 = tmp;
				signR = -signR;
			}

			for (j = 0, i = 0; i < fatherEdges_R; i ++)
			{
				if (j >= sonEdges) {
					cout << LINE_INFO << "Edges exceeds limit!" << endl;
					exit(0);
				}

				pCurEdgeInfo = &fatherVertex->GetGraph_R()[i];
				if (pCurEdgeInfo->edgenum == Ce)
				{
					if (GetPermNum(fatherVertex->GetGraph_R(), fatherEdges_R, Cn2))
						signR = -signR;
					continue;
				}

				// Enforce the effect of merging node (Cn2) to node (Cn1) and retaining (n1).
				sonEdgelist_R[j].node1 = (pCurEdgeInfo->node1 == Cn2) ? Cn1 : pCurEdgeInfo->node1;
				sonEdgelist_R[j].node2 = (pCurEdgeInfo->node2 == Cn2) ? Cn1 : pCurEdgeInfo->node2;

				sonEdgelist_R[j].edgenum = pCurEdgeInfo->edgenum;
				sonEdgelist_R[j].isSrc = pCurEdgeInfo->isSrc;
				sonEdgelist_R[j].mark = 0;
				j ++;

			}
			sonNodes = fatherNodes - 1;
			return sonEdgelist_R;
		}
		else if (Stype == CS && Ctype == VC)  // (G) : open CS, open VC
		{
#ifdef  TRACE			
			cout << LINE_INFO << "\n\tWorking on a VC-CS edge for [SYMBOL_OUT]" << endl;
#endif			

			for (j = 0, i = 0; i < fatherEdges_R; i ++)
			{
				if (j >= sonEdges) {
					cout << LINE_INFO << "Edges exceeds limit!" << endl;
					exit(0);
				}

				pCurEdgeInfo = &fatherVertex->GetGraph_R()[i];
				if (pCurEdgeInfo->edgenum == Ce)
					continue;
				
				sonEdgelist_R[j].node1 = pCurEdgeInfo->node1;
				sonEdgelist_R[j].node2 = pCurEdgeInfo->node2;
				
				sonEdgelist_R[j].edgenum = pCurEdgeInfo->edgenum;
				sonEdgelist_R[j].isSrc = pCurEdgeInfo->isSrc;
				sonEdgelist_R[j].mark = 0;
				j ++;
			}
			// Number nodes unchanged.
			return sonEdgelist_R;
		}
		else if (Stype == VS && Ctype == CC)  // (H) : short VS, short CC
		{
#ifdef  TRACE			
			cout << LINE_INFO << "\n\tWorking on a CC-VS edge for [SYMBOL_OUT]" << endl;
#endif			

			if (Cn1 == Sn1 && Cn2 == Sn2)
			{
				delete[] sonEdgelist_R;
				sonEdgelist_R = 0;
				return NULL;
			}

			// The controlling edge is to be shorted
			if (Cn1 > Cn2)
			{
				tmp = Cn2,    Cn2 = Cn1,    Cn1 = tmp;
				signR = -signR;
			}
			
			// The source edge is to be shorted
			if (Sn1 > Sn2)
			{
				tmp = Sn2,    Sn2 = Sn1,    Sn1 = tmp;
				signR = -signR;
			}

			if (Sn2 == Cn2)
			{
				if (Sn1 < Cn1)
					n1 = Sn1,    n2 = Sn1;
				else
					n1 = Cn1,    n2 = Cn1;	// (* corrected *)
					//n2 = Cn1;	//why not n2 = Cn1?	
			}
			else if (Sn2 == Cn1)
				n1 = Sn1,    n2 = Sn1;
			else if (Sn1 == Cn2)
				n1 = Cn1,    n2 = Cn1;
			else
				n1 = Sn1,    n2 = Cn1;

			for (j = 0, i = 0; i < fatherEdges_R; i ++)
			{
				if (j >= sonEdges) {
					cout << LINE_INFO << "Edges exceeds limit!" << endl;
					exit(0);
				}

				pCurEdgeInfo = &fatherVertex->GetGraph_R()[i];
				if (pCurEdgeInfo->edgenum == Se)
				{
					if (GetPermNum(fatherVertex->GetGraph_R(), fatherEdges_R, Sn2))
						signR = -signR;
					continue;
				}
				else if (pCurEdgeInfo->edgenum == Ce)
				{
					if (GetPermNum(fatherVertex->GetGraph_R(), fatherEdges_R, Cn2))
						signR = -signR;
					continue;
				}

				// Enforce the effect of merging node (Sn2) to node (Sn1) and retaining (n1).
				sonEdgelist_R[j].node1 = (pCurEdgeInfo->node1 == Sn1 || pCurEdgeInfo->node1 == Sn2) ? 
									n1 : pCurEdgeInfo->node1;
				sonEdgelist_R[j].node2 = (pCurEdgeInfo->node2 == Sn1 || pCurEdgeInfo->node2 == Sn2) ? 
									n1 : pCurEdgeInfo->node2;

				sonEdgelist_R[j].edgenum = pCurEdgeInfo->edgenum;
				sonEdgelist_R[j].isSrc = pCurEdgeInfo->isSrc;
				sonEdgelist_R[j].mark = 0;
				j ++;

			}

			for (i = 0; i < sonEdges; i ++)
			{
				pCurEdgeInfo = &sonEdgelist_R[i];
				if (pCurEdgeInfo->node1 == Cn2 || pCurEdgeInfo->node1 == Cn1)
					pCurEdgeInfo->node1 = n2;
				
				if (pCurEdgeInfo->node2 == Cn2 || pCurEdgeInfo->node2 == Cn1)
					pCurEdgeInfo->node2 = n2;
			}
			sonNodes = fatherNodes - 2;
			return sonEdgelist_R;
		}
		else 
		{
			delete[] sonEdgelist_R;
			sonEdgelist_R = 0;
			return NULL;
		}
	}
}


int GPDD::CreateSymbList()
{
	cout << LINE_INFO << endl;
	
	// Need an array for edge marks
	int *edge_mark = new int[NumEdges];

	int edge = 0, e = 0, j = 0;		

	for (e = 0; e < NumEdges; e++)
		edge_mark[e] = 0;

	// Traverse the Edge List.

	// Skip the first two edges which are a pair of I/O.
	//for (e = 2; e < E; e++)
	//Change by Xu Hui
	for (e = 0; e < NumEdges; e++)
	{
		// Skip the marked edge.
		if ( edge_mark[e] )		
			continue;

		switch (EdgeList[e].type)
		{
		case  Y:
		case  Z:
			edge = e;
			break;
		case  CC:
		case  VC:
			edge = EdgeList[e].pnum;	// symbol with the pairing edge
			edge_mark[edge] = 1;	// mark the pairing edge
			DpenSrcN ++;
			break;
		case  CS:
		case  VS:
			edge = e;			// symbol with this edge
			edge_mark[ EdgeList[e].pnum ] = 1;	// mark the pairing edge
			DpenSrcN ++;	
			break;
		default:
			cout << LINE_INFO << "Error: no such kind of edge (" << EdgeList[e].name << endl;
			return (STOP);
		}
		
		AddSymbol( &EdgeList[edge]);
		// Add symbol to the SymbListTemp[ ] where "lumping" has been taken care of
		// by creating sublists. Hence, the SymbListTemp[ ] is already lumped.
		
	}
	
	// Add the remaining one symbols.
	// Add "+" and "-" to the high end of SymbList[].
	SymbNode *pNode = SymbListTemp;
	SymbNode *pNewNode;

	LumpedSymbList = new SymbNode*[NumSymbols];

	for (int i = 0; i < NumSymbols; i++)
	{
		LumpedSymbList[i] = pNode;
		LumpedSymbList[i]->symbindex = i;  // save the symbol index 

		//Assign sensitivity to the symbol
		//LumpedSymbList[i]->Sens = LumpedSymbList[i]->edge->Sens;
		pNode = pNode->next;
	}

	cout << LINE_INFO << "Sort the symbol list" << endl;
	SortSymbList();	//*** Doing nothing there ...

	// Count the total number of symbols ...	
	TotalNumSymbols = 0;
	pNode = SymbListTemp;
	SymbNode *pParaNode;

	for (int i = 0; i < NumSymbols; i++)
	{
		TotalNumSymbols ++;
		
		pParaNode = LumpedSymbList[i]->nextP;

		if (pParaNode)
		{
			// The sublist exists. Then skip the first element in the substring,
			// which is duplicated by the ^element.
			pParaNode = pParaNode->nextP;
		}

		for ( ; pParaNode; pParaNode = pParaNode->nextP) {
			TotalNumSymbols ++;
		}  
	}

	/** The difference between ExpSymbList[] and FullSymbList[]:
		In SExapList[], the listing symbols (including the ^elements) are
			duplicated, then the sublist elements (skipping the ones already 
			included under the name of ^element) are appended to the expanded list.
		In FullSymbList[], the listing symbols (with the ^elements replaced by 
			its real elements) are duplicated, then the sublist elements 
			(skipping the ones already included) are appended to the expanded list.
		(G. Shi, 09-02-2010)
	*/

	ExpSymbList = new SymbNode*[TotalNumSymbols];
	
	// Scan the Lumped symbol list to process the listed elements.
	for (int i = 0; i < NumSymbols; i++)
	{
		ExpSymbList[i] = LumpedSymbList[i];	// the ^elements are duplicated as well
		ExpSymbList[i]->symbindex0 = i;
	}

	// Rescan the Lumped symbol list to process the sublists and do appending.
	j = NumSymbols;
	for (int i = 0; i < NumSymbols; i++)
	{
		pParaNode = LumpedSymbList[i]->nextP;

		if (pParaNode)
		{
			pParaNode = pParaNode->nextP;	// skip the sublist leader
		}

		while (pParaNode)
		{
			ExpSymbList[j] = pParaNode;
			ExpSymbList[j]->symbindex0	= j;
			ExpSymbList[j]->isSrc  = 0;
			ExpSymbList[j]->pairedge  = NULL;
			pParaNode = pParaNode->nextP;
			j ++;
		}
	}

/*	for(i = 0; i < TotalNumSymbols; i ++)
	{
		cout << "symbol " << i << "  " << ExpSymbList[i]->symbname << " " << SExendl;
	}*/

//	delete []LumpedSymbList;

	cout << LINE_INFO << "Total number of symbols : " << j << endl;

	cout << LINE_INFO << "Print the symbol list" << endl;
	PrintSymbols();

	LumpedSymbList = ExpSymbList;

	FullSymbList = new SymbNode*[TotalNumSymbols];
	
	// Scan the Lumped symbol list to process the listed leading elements.
	for (int i = 0; i < NumSymbols; i++)
	{
		// If there exists a sublist, use the sublist leader rather then the ^element.
		if (LumpedSymbList[i]->symbname[0] != '^')
			FullSymbList[i] = LumpedSymbList[i];
		else
		{
			FullSymbList[i] = LumpedSymbList[i]->nextP;
			FullSymbList[i]->isSrc = LumpedSymbList[i]->isSrc;
		      FullSymbList[i]->pairedge = LumpedSymbList[i]->pairedge;
		}
	}

	// Rescan the Lumped symbol list to process the sublists and do appending.
	j = NumSymbols;
	for (int i = 0; i < NumSymbols; i++)
	{
		pParaNode = LumpedSymbList[i]->nextP;

		if (pParaNode)
		{
			pParaNode = pParaNode->nextP;	// skip the sublist leader
		}

		while (pParaNode)
		{
			FullSymbList[j] = pParaNode;
			FullSymbList[j]->symbindex0	= j;
			FullSymbList[j]->isSrc = 0;
			FullSymbList[j]->pairedge = NULL;
			pParaNode = pParaNode->nextP;
			j ++;
		}
	}

	//cout<<LINE_INFO<<"Print out No Lump List"<<endl;
	for (int i = 0; i < TotalNumSymbols; i ++)
	{
		/* cout << FullSymbList[i]->symbname << " " << FullSymbList[i]->edge->node1 << " " 
			<< FullSymbList[i]->edge->node2 <<endl; */
	
		// Added by Xu Hui
		// Create symbol name to SymbNode hash table.
		//SymbHashTable[string(FullSymbList[i]->symbname)] = &(FullSymbList[i]->value);
		
		SymbHashTable[string(FullSymbList[i]->symbname)] = FullSymbList[i];
	}

	delete [] edge_mark;
	edge_mark = 0;
	
	return (OK);
}


int GPDD::DetectTree(EdgeInfo *pEdgeInfo, int e_size, int n_size)
{
	int 	nodes = 0;

	if (e_size == 0)
		return 1;

	short *mark = new short[NumNodes];

	memset(mark, 0, NumNodes* sizeof(short));
	
	for (int i = 0; i < e_size; i ++)
	{
		if (!mark[pEdgeInfo[i].node1])
		{
			mark[pEdgeInfo[i].node1] = 1;
			nodes ++;
		}
		if (!mark[pEdgeInfo[i].node2])
		{
			mark[pEdgeInfo[i].node2] = 1;
			nodes ++;
		}
	}
	
	delete [] mark;
	mark = 0;

	if (nodes == n_size)
		return 1;

	return 0;
}


/**
  Get the symbol value given its name.
  Added by G. Shi.
*/
int GPDD::GetSymbolValue(const string &symbname, double &value)
{
	
	if (NumSymbols == 0)
	{
		cout << LINE_INFO << "\nEmpty symbol list.\n" << endl;
		return (NOPE); // an error value
	}
	
	for (int n=0; n < NumSymbols; ++n)
	{
		SymbNode *pNode = LumpedSymbList[n];

		for ( ; pNode; pNode = pNode->nextP)
		{
			string name(pNode->symbname);
			if (name == symbname)
			{
				value = pNode->value;
				cout << LINE_INFO << "\n\t" 
					<< name << " : value (" << pNode->value << ") " << endl;
				return (OK);
			}
		}
	}
	
	cout << LINE_INFO << "Symbol \" " << symbname << "\" not found." << endl;
	return (NOPE);
}

/**
	Get the symbol value given its index and name. 
	(Function name changed by G.Shi, 08-30-2010).
*/
//double GPDD::Search_value(int index, char * name )
double GPDD::GetSymbolValue(int index, const char * name )
{
	double intial_value = 0.0;
	SymbNode *pNode;
	
	if(strcmp(LumpedSymbList[index]->symbname, name) == 0)
	{
		intial_value = LumpedSymbList[index]->value;
	}
	else
	{
		pNode = LumpedSymbList[index]->nextP;
		for(; pNode; pNode = pNode->nextP)
		{
			if(strcmp(pNode->symbname, name) == 0)
			{
				intial_value = pNode->value;
				break;
			}
		}
	}
	
	return intial_value;
}


#if 00
/**
	Set the value of the symbol under certain frequency.
	<BR> R : value = 1/R;
	<BR> G : value = G;
	<BR> C : value = 0 + C * 2 * pi * freq * j;
	<BR> L : value = 0 - 1/(L * 2 * pi * freq)j;
	<BR> sources: value = the gain of the sources;
	@param freq the frequency.
*/

void GPDD::EvaluateSymbolAtFreq(double freq)
{
    SymbNode	*pNode;
    double		real = 0, imag = 0;		
    double		value;

    //for sensitivity
    double sr=0, si=0;
    double sens;
	
    for (int i = 0; i < NumSymbols; i++)
    {
	if ((pNode = LumpedSymbList[i]->nextP))	// Has lumped edges.
	{
	    real = imag = 0;		// initialize
	    for (; pNode; pNode = pNode->nextP)
	    {
		value = pNode->value;

		//Add sensitivity if lumped
		sens = pNode->sens_r;

		switch (pNode->edge->name[0])
		{
		    case 'R':
			real += (value == 0? DBL_MAX : 1/value);
			sr += (value == 0? DBL_MAX : -sens/(value*value));
			break;
		    case 'G':
			real += value;
			sr += sens;
			break;
		    case 'L':
			//Why "-="? Not 1/Ls but -1/Ls?
			imag -= (value == 0? DBL_MAX : 1/(2*PI*freq*value) );
			si += (value == 0? DBL_MAX : -sens/(2*PI*freq*value*value));
			break;
		    case 'C':
			imag += (2*PI*freq*value);
			si += sens*2*PI*freq;
			break;
		    default:
			break;
		}
	    }
	    // Store the real and imag values.
	    LumpedSymbList[i]->value = real;
	    LumpedSymbList[i]->imag  = imag;

	    // Store the sensitivity
	    LumpedSymbList[i]->sens_r = sr;
	    LumpedSymbList[i]->sens_i = si;
	}
    }
}
#endif


/**
	Recalculate all symbol values given a new frequency.
*/
void GPDD::EvaluateAllSymbolsAtFreq(double freq)
{
	SymbNode	*pNode;
	double		real = 0, imag = 0;		
	double		value;

	for (int i = 0; i < NumSymbols; i++)
	{
		if ((pNode = LumpedSymbList[i]->nextP))	// Has lumped edges.
		{
			real = imag = 0;		// initialize; used for accumulating the sublist values
			
			for (; pNode; pNode = pNode->nextP)
			{
				value = pNode->value;

				switch (pNode->edge->name[0])
				{
				case 'R':
					real += (value == 0? DBL_MAX : 1/value);
					break;
				case 'G':
					real += value;
					break;
				case 'L':
					imag -= (value == 0? DBL_MAX : 1/(2*PI*freq*value) );
					break;
				case 'C':
					imag += (2*PI*freq*value);
					break;
				default:
					break;
				}
			}
			
			// Store the real and imag values.
			LumpedSymbList[i]->value = real;
			LumpedSymbList[i]->imag  = imag;
		}
	}
}

/**
	Expand the composite GPDD nodes (corresponding to the lumped symbols)
	to individual GPDD nodes.

	This function must be called before calling sExpandGPDD()
	to derive the s-coefficients.
*/
//void GPDD::sexpand()	// Old name by W. Chen
int  GPDD::ExpandLumpedGPDD()
{
	cout << LINE_INFO << endl;

	if (LumpedExpandedFlag) {
		cout << LINE_INFO << "\n\tThe GPDD has been lump-expanded." << endl;
		return (OK);
	}
	

///	UnMarkAllNodes(pGPDDroot, 1);

	if (!pGPDDroot) {
		cout << LINE_INFO << "GPDD root is NULL." << endl;
		return (STOP);
	}

	ExpOne = new ExpGPDDnode;
	ExpZero = new ExpGPDDnode;

	ExpOne->pGPDDnode = VertexOne;
	ExpOne->SymbolIndex = -1;
	ExpOne->Sign_OUT     = 1;
	ExpOne->Sign_IN     = 1;
	ExpOne->Vertex_OUT     = NULL;
	ExpOne->Vertex_IN    = NULL;

	ExpZero->pGPDDnode = VertexZero;
	ExpZero->SymbolIndex = -2;
	ExpZero->Sign_OUT     = 1;
	ExpZero->Sign_IN     = 1;
	ExpZero->Vertex_OUT     = NULL;
	ExpZero->Vertex_IN    = NULL;

	TotalExpandedGPDDnodes  = 0;

	TraversalMark = pGPDDroot->GetMark();
	FlipTraversalMark();
	
	ExpGPDDroot = ExpandLumpedGPDD_R(pGPDDroot);
	
	pGPDDroot->SetMark(TraversalMark);

	LumpedExpandedFlag = true;	// have to expand the s-polynomial only once

	return (OK);
}

	
ExpGPDDnode* GPDD::ExpandLumpedGPDD_R(GPDDnode* p)
{
	if (p->GetMark() == TraversalMark)
	{
		return p->pExpandNode;
	}
	else
	{
		p->SetMark(TraversalMark);
	}

	int tag;
	ExpGPDDnode *pNewNode = new ExpGPDDnode();
	ExpGPDDnode *pSharedNode;

	pNewNode->pGPDDnode = p;
	pNewNode->SymbolIndex = p->GetSymbolIndex();
	pNewNode->Sign_OUT = p->GetSign_OUT();
	pNewNode->Sign_IN = p->GetSign_IN();

	pSharedNode = Share_ExpandedGPDDnode(pNewNode, tag);

	if (tag)
	{
		p->pExpandNode = pSharedNode;
		return pSharedNode;
	}

	ExpGPDDnode *pShort, *pOpen;

	if (p->GetVertex_IN() == VertexOne)
	{
		pShort = ExpOne;
	}
	else if (p->GetVertex_IN() == VertexZero)
	{
		pShort = ExpZero;
	}
	else
	{
		pShort = ExpandLumpedGPDD_R(p->GetVertex_IN());
	}


	if (p->GetVertex_OUT() == VertexOne)
	{
		pOpen = ExpOne;
	}
	else if (p->GetVertex_OUT() == VertexZero)
	{
		pOpen = ExpZero;
	}
	else
	{
		pOpen  = ExpandLumpedGPDD_R(p->GetVertex_OUT());
	}

	int 	signS, signO;
	SymbNode *pCurSymb;
	signS = p->GetSign_IN();
	signO = p->GetSign_OUT();

	if (LumpedSymbList[p->GetSymbolIndex()]->symbname[0] == '^')
	{
		pSharedNode->Sign_OUT = 1;
	//	totalExpandedNodes ++ ;
	//	cout << totalExpandedNodes << " " << LumpedSymbList[p->GetSymbolIndex()]->symbname << endl;
		pSharedNode->Vertex_IN = pShort;
		
		pCurSymb = LumpedSymbList[p->GetSymbolIndex()];

		pCurSymb = pCurSymb->nextP->nextP;

		pNewNode = new ExpGPDDnode(pShort, pOpen, NULL, 
			pCurSymb->symbindex0, p->GetSign_IN(), p->GetSign_OUT());

		pSharedNode->Vertex_IN = pShort;
		pSharedNode->Vertex_OUT  = pNewNode;

		pCurSymb = pCurSymb->nextP;

		while (pCurSymb)
		{
			pNewNode = new ExpGPDDnode(pShort, pSharedNode->Vertex_OUT, NULL, 
				pCurSymb->symbindex0, p->GetSign_IN(), 1);
			TotalExpandedGPDDnodes ++;
			pSharedNode->Vertex_OUT = pNewNode;
			pCurSymb = pCurSymb->nextP;
		}

		p->pExpandNode = pSharedNode;
		return pSharedNode;

	}

	pSharedNode->Vertex_IN = pShort;
	pSharedNode->Vertex_OUT  = pOpen;

	p->pExpandNode = pSharedNode;
	return pSharedNode;
	
}


EdgeInfo* GPDD::InitLeftGraph(int &signL)
{
	cout << LINE_INFO << endl;

	if (!EdgeList) {
		cout << LINE_INFO << "Error : No edge list." << endl;
		exit(0);
	}

	EdgeInfo *edgelist_L;

	NumEdges_L = NumSymbols + DpenSrcN - VcN;

	if (NumEdges_L <= 0) {
		cout << LINE_INFO << "\n\tIncorrect number of edges in the Left Graph."
			<< "\n\tProgram is terminated." << endl;
		exit(0);
	}
	
	edgelist_L = new EdgeInfo[NumEdges_L];

	int	j = 0;
	for (int i = 0; i < NumEdges; i ++)
	{
		if (EdgeList[i].lumped || EdgeList[i].type == VC)
			continue;

		edgelist_L[j].node1 = EdgeList[i].node1;
		edgelist_L[j].node2 = EdgeList[i].node2;
		edgelist_L[j].mark = 0;
		edgelist_L[j].edgenum = i;
		edgelist_L[j].isSrc = EdgeList[i].type > 1 ? 1 : 0;
		j ++;
	}

	cout << LINE_INFO << "NumEdges_L = " << NumEdges_L 
		<< ", j = " << j << endl;

	signL = 1;
	
	return edgelist_L;
}

	
EdgeInfo* GPDD::InitRightGraph(int &signR)
{
	cout << LINE_INFO << endl;

	EdgeInfo* edgelist_R;

	NumEdges_R = NumSymbols + DpenSrcN - CsN;

	if (NumEdges_R <= 0) {
		cout << LINE_INFO << "\n\tIncorrect number of edges in the Right Graph."
			<< "\n\tProgram is terminated." << endl;
		exit(0);
	}

	edgelist_R = new EdgeInfo[NumEdges_R];

	int 	j = 0;
	for (int i = 0; i < NumEdges; i ++)
	{
		if (EdgeList[i].lumped || EdgeList[i].type == CS)
			continue;

		edgelist_R[j].node1 = EdgeList[i].node1;
		edgelist_R[j].node2 = EdgeList[i].node2;
		edgelist_R[j].mark = 0;
		edgelist_R[j].edgenum = i;
		edgelist_R[j].isSrc = EdgeList[i].type > 1 ? 1 : 0;
		j ++;
	}

	cout << LINE_INFO << "NumEdges_R = " << NumEdges_R 
		<< ", j = " << j << endl;
	
	signR = 1;

	return edgelist_R;
}

/**
  Initialize the GPDD root vertex.
*/
void GPDD::InitGPDDroot()
{
	cout << LINE_INFO << endl << flush;
	
	EdgeInfo 		*edgelist_L, *edgelist_R;
	Subgraph 	*graph_L, *graph_R, *hashgraph_L, *hashgraph_R;

	int signL = 1; 
	int signR = 1;
	int hashed;
	
	edgelist_L = InitLeftGraph(signL);
	edgelist_R = InitRightGraph(signR);
	// NumEdges_L & NumEdges_R are fixed in the above calls

//	int	numedges_L = NumSymbols + DpenSrcN - VcN;
//	int	numedges_R = NumSymbols + DpenSrcN - CsN;

	cout << LINE_INFO << "\n\tCheck the initial graph-pair --" << endl;
	cout << LINE_INFO << "Edges of the (initial) L-graph: " << endl;
	PrintPartialEdgeList(edgelist_L, NumEdges_L);

	cout << LINE_INFO << "Edges of the (initial) R-graph: " << endl;
	PrintPartialEdgeList(edgelist_R, NumEdges_R);


	//Add NumNodes to Subgraph
	graph_L = new Subgraph(edgelist_L, NumEdges_L, NumNodes, signL, NumNodes);
	graph_R = new Subgraph(edgelist_R, NumEdges_R, NumNodes, signR, NumNodes);

	hashgraph_L = ShareSubgraph(graph_L);
	hashgraph_R = ShareSubgraph(graph_R);
	
	pGPDDroot = new GPDDnode(LumpedSymbList[0], 1, 1, hashgraph_L, hashgraph_R); 

	ShareGPDDnode(pGPDDroot, hashed);
	
}

int GPDD::GetNodeSize(EdgeInfo * pEdgeInfo, int e_size)
{
	int i, n1, n2, totalnodes;
	memset(NodeMark, 0, NumNodes* sizeof(int));
	
	for(i = 0; i < e_size; i ++)
	{
		n1 = pEdgeInfo[i].node1;
		n2 = pEdgeInfo[i].node2;
		NodeMark[n1] = 1;
		NodeMark[n2] = 1;
	}
	
	totalnodes = 0;

	for(i = 0; i < NumNodes; i ++)
	{
		if(NodeMark[i])
			totalnodes ++;
	}

	return totalnodes;
}

/**
	NR: non-recursive implementation of graph reduction. 
*/
void GPDD::GraphReduction_NR()
{
	cout << LINE_INFO << endl;
	
	GPDDnode *pCurNode, *pNewNode, *pHashNode;

	stack<GPDDnode*> *NodeStack = new stack<GPDDnode*>;

	if (!NodeStack) {
		cout << LINE_INFO << "Stack allocation error." << endl;
		exit(0);
	}

	TotalGPDDnodes = 0;
	GPDDnodeSharingCount = 0;
	int hashed = 0;
	
	if (!pGPDDroot)
		return;
	
	NodeStack->push(pGPDDroot);

	while (!NodeStack->empty())
	{
//		cout << LINE_INFO << "Code is here." << endl << flush;

		if (NodeCreationCount > 0 && NodeCreationCount % 10000 == 0)
		{
			cout << LINE_INFO << "NodeCreationCount = " << NodeCreationCount << endl; 
		} 
		
		pCurNode = NodeStack->top();
		
//		cout << LINE_INFO << "Stack top is [" << pCurNode->GetSymbol()->symbname 
//			<< "]" << endl; 

		if (!pCurNode->GetVertex_IN())
		{
//			cout << LINE_INFO << "Code is here." << endl << flush;

			NodeCreationCount ++;
			pNewNode = CreateGPDDnode(pCurNode, SYMBOL_IN);

	
			if (pNewNode == VertexZero || pNewNode == VertexOne || !pNewNode)
			{
				pCurNode->SetVertex_IN(pNewNode);
				//cout << LINE_INFO << "Continue without hashing ---" << endl;
				continue;
			}

//			cout << LINE_INFO << "Code is here." << endl << flush;
			
			pHashNode = ShareGPDDnode(pNewNode, hashed);
			if (!hashed) {
				NodeStack->push(pHashNode);	// This is a newly created <IN> node.
			}
			
			pCurNode->SetVertex_IN(pHashNode);
			
		}
		else if (!pCurNode->GetVertex_OUT())
		{
//			cout << LINE_INFO << "Code is here." << endl << flush;

			NodeCreationCount ++;
			pNewNode = CreateGPDDnode(pCurNode, SYMBOL_OUT);
	
			if (pNewNode == VertexZero || pNewNode == VertexOne || !pNewNode)
			{
				pCurNode->SetVertex_OUT(pNewNode);
				//cout << LINE_INFO << "Continue without hashing ---" << endl;
				continue;
			}
			
			pHashNode = ShareGPDDnode(pNewNode, hashed);
			if (!hashed)
			{
				NodeStack->push(pHashNode);	// This is a newly created <OUT> node.
			}
			
			pCurNode->SetVertex_OUT(pHashNode);
			
		}
		else
		{
//			cout << LINE_INFO << "Code is here." << endl << flush;

			NodeStack->pop();	// Node with both children created is popped and counted.
			TotalGPDDnodes ++;
		}

#ifdef	TRACE		
		if (TotalGPDDnodes % 10000 == 0)
		{
			cout << LINE_INFO << "TotalGPDDnodes = " << TotalGPDDnodes << endl; 
		} 
#endif
	}

	delete NodeStack;
	NodeStack = 0;

	SubgraphCount = 0;

	cout << LINE_INFO << "\n\tDeleting the uniqSubgraph hash table ..." << endl;
	
	if (uniqSubgraph)
	{
		for (int i = 0 ; i < SUBGRAPH_PRIME ; i ++ )
		{
			if (!uniqSubgraph[i])
				continue;
			
			SubgraphCount += uniqSubgraph[i]->Size();
			
			uniqSubgraph[i]->SetDeleteFlag(TRUE, FALSE);    
			delete uniqSubgraph[i];
			uniqSubgraph[i] = NULL;
		}
		
		delete []uniqSubgraph;
		uniqSubgraph = NULL;
	}
}

/**
  The key member function of GPDD, which calls the recursive implementations.
*/
void GPDD::GraphReduction()
{
	cout << LINE_INFO << endl;
	
	GPDDnode *vertex_out, *vertex_in;
	
	if (NumNodes > NumEdges + 1)
//	if (NumNodes > NumEdges + 2)		//Since GND has been added to the graph in HybridSim (Xu Hui)
	{
		// Graph is disconnected under this condition
		pGPDDroot = NULL;
		return;
	}
	
	NodeCreationCount = 0;
	SubgraphSharingCount = 0;
	GPDDnodeSharingCount = 0;

	//InitGPDDroot();
	// It is now initialized by NPortGPDD.

	// Enter recursive Binary Decisions
	cout << LINE_INFO << "Start recursion on <INCLUDE> ... " << endl;
	vertex_in = Reduction_INCLUDE_R(pGPDDroot);

	cout << LINE_INFO << "Start recursion on <EXCLUDE> ... " << endl;
	vertex_out = Reduction_EXCLUDE_R(pGPDDroot);

	pGPDDroot->SetVertex_IN(vertex_in);
	pGPDDroot->SetVertex_OUT(vertex_out);

	cout << LINE_INFO << "-- DONE" << endl;
}

void GPDD::PrintGPDD(ostream &out)
{
	cout << LINE_INFO << endl;
	
//	UnMarkAllNodes(pGPDDroot, 1);

	cout << "\n  Notation:" << endl;
	cout << "\t[#index];  (left) / for VERTEX_IN; (right) \\ for VERTEX_OUT" << endl;
	cout << "\tat (address)" << endl;

	TraversalMark = pGPDDroot->GetMark();
	FlipTraversalMark();
	
	PrintGPDD_R(out, pGPDDroot);

	pGPDDroot->SetMark(TraversalMark);
}
	
void GPDD::PrintGPDD_R(ostream &out, GPDDnode *p)
{
	if (p->GetMark() == TraversalMark)
		return;
	
	p->SetMark(TraversalMark);

	p->PrintGPDDnode(out);

	if(p->GetVertex_IN())
		PrintGPDD_R(out, p->GetVertex_IN());
	if(p->GetVertex_OUT())
		PrintGPDD_R(out, p->GetVertex_OUT());
}

/**
	Print all the terms in the GRDD.
	Note: Memory for PathStack is created and deleted inside this function.
	(Same code as GRASS3, G. Shi, 09-03-2010)
*/
void GPDD::PrintTerms(ostream &out)
{
	if (NumSymbols == 0)
		return;
	
	if (pGPDDroot->GetSymbolIndex() < 0) {
		return;
	}
	
	if (PathStack) 
	{
		delete PathStack;		// Just to prevent memory leakage!
	}
	
	PathStack = new int[NumSymbols];
	if (!PathStack) {
		cout << LINE_INFO << "Memory Allocation Failure!" << endl;
		exit(0);
	}
	
	NumTerms = 0;
	TermSign	= 1;
	StackTop = 0;

	printf("\n[%s : line %d] %s() : Printing GRDD terms ... \n", __FILE__, __LINE__, __FUNCTION__);
	
	PrintTerms_R(out, pGPDDroot);

	delete [] PathStack;
	PathStack = NULL;
}


/*
    The recursion function for PrintPath().
    We traverse the GRDD from the root. Everytime seeing the vertex One, 
    we print the term. An array is used and the symbol indexes along the one-path are
    stored in the array for priniting the symbol names.
	(Same code as GRASS3, G. Shi, 09-03-2010)
*/

void GPDD::PrintTerms_R(ostream &out, GPDDnode *p)
{
	int 	sign = 1;

	if (p->GetSymbolIndex() == -1)	// Reached vertex [One]
	{
		out << "Term: \t";	
		if (TermSign < 0)
			out << "-";
		else
			out << "+";

		for (int i = 0; i < StackTop; i ++)
		{
			if (LumpedSymbList[PathStack[i]]->edge->type == Z)
			{ 
				out << "(";	// erap the Z-type symbname by parentheses
				out << LumpedSymbList[PathStack[i]]->symbname;
				out << ").";
			}
			else
				out << LumpedSymbList[PathStack[i]]->symbname << ".";
		}
		out << endl;	

		NumTerms ++;
		return;
	}
	else if (p->GetSymbolIndex() == -2)	{
		// We only have two terminating vertices: One (-1) and Zero (-2)
		return;	// Reached one of the treminal vertices: One or Zero
	}
	
	sign = TermSign;	// Save the current sign
	
	PathStack[StackTop++] = p->GetSymbolIndex();	// save the current symbindex

	TermSign *= p->GetSign_IN();	// update the path sign

/*------------------------------------------------------------------*/	
	PrintTerms_R(out, p->GetVertex_IN());		// <-- [SYMBOL_IN] arrow
/*------------------------------------------------------------------*/	

	TermSign = sign * p->GetSign_OUT();
	StackTop --;		// Roll back the top pointer
	
/*------------------------------------------------------------------*/	
	PrintTerms_R(out, p->GetVertex_OUT());	// <-- [SYMBOL_OUT] arrow
/*------------------------------------------------------------------*/	
	
	TermSign = sign;	// Recover the TermSign
}



int GPDD::GetPermNum(EdgeInfo *pEdgeInfo, int e_size, int n)
{
	memset(PermMark, 0, NumNodes* sizeof(int));


	int i, permN = 0;

	for(i = 0; i < e_size; i ++)
	{
		if(pEdgeInfo[i].node1 < n && !PermMark[pEdgeInfo[i].node1])
		{
			PermMark[pEdgeInfo[i].node1] = 1;
			permN ++;
		}
		if(pEdgeInfo[i].node2 < n && !PermMark[pEdgeInfo[i].node2])
		{
			PermMark[pEdgeInfo[i].node2] = 1;
			permN ++;
		}
	}

	return (permN & 0x1);
}

int GPDD::CountValidEdges(EdgeInfo * pTree, int En)
{
	int	validedges = 0;
	
	for (int i = 0; i < En;  i ++)
	{
		if(pTree[i].node1 ^ pTree[i].node2)
			validedges ++;
	}

	return validedges;
}

/**
	GPDD reduce.
	After reducion, there will not be isomorphic subgraphs in GPDD.
	Isomorphic subgraphs is generated due to the isomorphic reduced graphs.
*/
int GPDD::Reduce()
{
	cout << LINE_INFO << endl;
	
	int i, j;
	GPDDnode_tag *ptag;
	GPDDnode *ps, *po;

	HashTree ***hashtree = new HashTree **[NumSymbols];

	NumNodesRemoved  = 0;

	// Initialize the hashtree table.
	for(i = 0; i < NumSymbols; i ++)
	{
		hashtree[i] = new HashTree*[TAG_PRIME2];
		for(j = 0; j < TAG_PRIME2; j ++)
			hashtree[i][j] = NULL;
	}

	VertexZero->GetTag()->SetTriple(1, 0, 0);
	VertexOne->GetTag()->SetTriple(2, 0, 0);
	
	NodeList.push_back(VertexZero);
	NodeList.push_back(VertexOne);

	NumNodesRemaining = NodeList.size();

	CollectNodes();			// Build the 2D ReduceList[]
	
	SetAllNodeLabels(hashtree);	// Build a new hash table with the help of GPDDnode_tags.

	if (pGPDDroot != NodeList[NumNodesRemaining - 1])
	{
		cout << LINE_INFO << "Reduce Error!" << endl;
		return (STOP);
	}

	/**
	  Construct the reduced GPDD using the GPDDnodes saved in the vector NodeList.
	  The parent-children relation has been saved in the tags, which are the indices
	  of NodeList[ ]. (See the function SetNodeLabel(). 
	*/
	for (i = 2; i < NumNodesRemaining; i ++)
	{
		ptag = NodeList[i]->GetTag();

		ps = NodeList[abs(ptag->label_left) - 1];
		po = NodeList[abs(ptag->label_right) - 1];
		
		NodeList[i]->SetVertex_OUT(po);
		NodeList[i]->SetVertex_IN(ps);
	}
	
	NodeList.clear();

	// Remove the hashtree
	for (i = 0; i < NumSymbols; i ++)
	{
		if (hashtree[i])
		{
			for (j = 0; j < TAG_PRIME2; j ++)
			{
				if (hashtree[i][j])
				{
					hashtree[i][j]->SetDeleteFlag(FALSE, FALSE);
					delete hashtree[i][j];
					hashtree[i][j] = 0;
				}
			}
			delete[] hashtree[i];
			hashtree[i] = 0;
		}
	}
	delete[] hashtree;
	hashtree = 0;

	return (OK);
}

void GPDD::GC()
{
	cout << LINE_INFO << "Running garbage collection ... " << endl;
	
	int i, j;
	for(i = GPDD_PRIME_LEVEL - 1; i >= 0; i --)
	{
		for(j = GPDD_PRIME - 1; j >= 0; j --)
		{
			if(!uniqGPDD[i][j])
				continue;
			if(uniqGPDD[i][j]->Size() != 0)  // unique[] is a list of HashTrees
	      			doGC(i, j);
		}
	}
}

void GPDD::doGC(int i, int j)
{
	GPDDnode *p;

	HashTree *unique_save = new HashTree();

	HashTree_For_Each_AnyObject(*uniqGPDD[i][j], key, object) {
		p = (GPDDnode *) key;
		if (p->NumReferences()) {
			unique_save->Put(p, *p);	// *** Put ***
		} else {
			ASSERT(!p->NumReferences(), "deleting a referenced node");
			delete p;	// Node released
			// This node is garbage collected.
		}
	} HashTree_For_Each_End;
	
	delete uniqGPDD[i][j];
	
	uniqGPDD[i][j] = unique_save;
}



void GPDD::PrintTransferFunction(ostream &out)
{
	cout << LINE_INFO << endl;
	
	if (NumSymbols == 0)
		return;
	
	if(pGPDDroot->GetSymbolIndex() < 0)
		return;
	
	if (PathStack) 
	{
		delete PathStack;		// Just to prevent memory leakage!
	}
	
	PathStack = new int[NumSymbols];
	if (!PathStack) {
		cout << LINE_INFO << "Memory Allocation Failure!" << endl;
		exit(0);
	}

	StackTop = 0;

	if (pGPDDroot->GetVertex_IN())
	{
		TermSign	= -1;
		
		cout << "\n  The Numerator of H(s):" << endl;
		StackTop = 0;
		PrintTransferFunction_R(out, pGPDDroot->GetVertex_IN());
	}

	if (pGPDDroot->GetVertex_OUT())
	{
		TermSign	= 1;
		
		cout << "\n  The Denominator of H(s):" << endl;
		StackTop = 0;
		PrintTransferFunction_R(out, pGPDDroot->GetVertex_OUT());
	}

	delete[] PathStack;
	PathStack = NULL;
}


	
void GPDD::PrintTransferFunction_R(ostream &out, GPDDnode *p)
{
	int sign = 1;
	
	if (p->GetSymbolIndex() == -1)
	{
		if (TermSign < 0)
			out << "    -";
		else
			out << "    +";

		for (int i = 0; i < StackTop; i ++)
		{
			if (LumpedSymbList[PathStack[i]]->edge->type == Z)
			{ 
				out << "(";
				out << LumpedSymbList[PathStack[i]]->symbname;
				out << ").";
			}
			else
				out << LumpedSymbList[PathStack[i]]->symbname << ".";
		}
		out << endl;	
		return;
	}
	
	if (p->GetSymbolIndex()< 0)
		return;
		
	sign = TermSign;
	PathStack[StackTop++] = p->GetSymbolIndex();

	TermSign *= p->GetSign_IN();	
	PrintTransferFunction_R(out, p->GetVertex_IN());

	TermSign = sign * p->GetSign_OUT();
	StackTop--;
	PrintTransferFunction_R(out, p->GetVertex_OUT());
	TermSign = sign;
}



/**
	Print all the coeffcients of the two s-polynomials of H(s) = N(s)/D(s).
	@param out the ostream object for console output.

	(Code from GRASS3, G. Shi, 09-03-2010)
*/
void GPDD::Print_sPolynomialCoeffs(ostream &out)
{
	cout << LINE_INFO << endl;
	
	int 		i;

	vector<sCoeffNode*> 	*pCoeff; 	// = ExpGPDDroot->Node_IN->sCoefList;
	vector<int> 			*pPower	; 	// = ExpGPDDroot->Node_IN->sPowerList;

	pCoeff  = ExpGPDDroot->Vertex_IN->sCoeffList;
	pPower = ExpGPDDroot->Vertex_IN->sPowerList;
	
	out << "\n  Numerator -- " << pPower->size() << " term(s)" << endl;

	for (i = 0; i < pCoeff->size(); i ++)
	{
		printf("\n    s^{%d}  ", (*pPower)[i]);
		PrintCoefficient(out, ((*pCoeff)[i]), ExpGPDDroot->Sign_IN);
	}

	pCoeff  = ExpGPDDroot->Vertex_OUT->sCoeffList;
	pPower = ExpGPDDroot->Vertex_OUT->sPowerList;
	
	out << "\n  Denominator -- " << pPower->size() << " term(s)" << endl;

	for (i = 0; i < pCoeff->size(); i ++)
	{
		printf("\n    s^{%d}  ", (*pPower)[i]);
		PrintCoefficient(out, ((*pCoeff)[i]), ExpGPDDroot->Sign_OUT);
	}
	
}



void GPDD::PrintCircuitStatistics()
{
#ifdef	TRACE
	printf("\n[%s : line %d] %s() ...", __FILE__, __LINE__, __FUNCTION__);
#endif

	printf("\n[%s : line %d] Graph statistics", __FILE__, __LINE__);
	cout << "\n\t# Edges    \t= " << GetNumEdges() 
		 << "\n\t# Nodes   \t= " << GetNumNodes()
		 << "\n\t# Symbols \t= " << GetNumSymbols()
		 <<"\n"
		 << "\n\t# VC edges = " << GetNumVC()
	       << "\n\t# VS edges = " << GetNumVS()
		 << "\n\t# CC edges = " << GetNumCC()
		 << "\n\t# CS edges = " << GetNumCS()
		 << "\n\n\t# Dep srcs = " << GetNumDependentSource() << "\n\n";
}

/**
	Print all the coeffcients of one s-exponent.
*/
void GPDD::PrintCoefficient(ostream &out, sCoeffNode *coefnode, int sign)
{
	if (NumSymbols == 0)
		return;
	
	if (coefnode->pExpNode->SymbolIndex < 0)
		return;
	
	if (PathStack) 
	{
		delete PathStack;		// Just to prevent memory leakage!
	}
	
	PathStack = new int[NumSymbols];
	if (!PathStack) {
		cout << LINE_INFO << "Memory Allocation Failure!" << endl;
		exit(0);
	}
	
	NumTerms = 0;
	TermSign	= sign;
	StackTop	= 0;

	PrintCoefficient_R(out, coefnode);

	delete[] PathStack;
	PathStack = NULL;
}



/**
	Recursively Print all the terms in the coeffcients of one s-exponent.
*/
void GPDD::PrintCoefficient_R(ostream &out, sCoeffNode *coefnode)
{
	int 		sign = 1, i;
	
	if (coefnode->pExpNode->SymbolIndex == -1)
	{
		out << "\n    Term:  ";
		
		if (TermSign < 0)
			out << "-";
		else
			out << "+";

		for (i = 0; i < StackTop; i ++)
		{
			if (LumpedSymbList[PathStack[i]]->edge->type == Z)
			{ 
				out << "(";
				out << LumpedSymbList[PathStack[i]]->symbname;
				out << ").";
			}
			else
				out << LumpedSymbList[PathStack[i]]->symbname << ".";
		}
		out << endl;	

		NumTerms ++;
		return;
	}
	else if (coefnode->pExpNode->SymbolIndex == -2)
		return;
	
	sign = TermSign;
	PathStack[StackTop++] = coefnode->pExpNode->SymbolIndex;

	TermSign *= coefnode->pExpNode->Sign_IN;	
	
	PrintCoefficient_R(out, coefnode->Vertex_IN);

	TermSign = sign * coefnode->pExpNode->Sign_OUT;
	StackTop --;
	
	PrintCoefficient_R(out, coefnode->Vertex_OUT);
	
	TermSign = sign;
}


/**	
	PrintFullEdgeList()  -- (newly added by G. Shi, 09-02-2010)
*/
void GPDD::PrintFullEdgeList()
{
	cout << LINE_INFO << endl;

	if (NumEdges == 0) {
		cout << LINE_INFO << "\n\tNo edge in the graph." << endl;
	}

	for (int i = 0; i < NumEdges; i ++)
	{
		printf("\n\t[%2d] \"%6s\"\t (%2d) ---- (%2d)  : paired with - %s", 
				i, 
				EdgeList[i].name,  	
				EdgeList[i].node1,
				EdgeList[i].node2, 
				(EdgeList[i].pname? EdgeList[i].pname : "(null)"));
	}
	cout << endl;
}

void GPDD::PrintGraph(ostream &out)
{
	out << LINE_INFO << "\n\tNumEdges = " << NumEdges 
		<< ", NumNodes = " << NumNodes << endl;
	
	out << setiosflags(ios::right)
		 << setw(10) << "edge#"
		 << setw(10) << "edge"
		 << setw( 5) << "n1#"
		 << setw( 5) << "n2#"
		 << setw( 6) << "type"
		 << setw(10) << "pair"
		 << setw(20) << "value"
		 << endl;

	for (int e = 0; e < NumEdges; e++)
	{
		out << setiosflags(ios::right)
			<< setw(10) << e
			<< setw(10) << EdgeList[e].name 
			<< setw( 5) << EdgeList[e].node1
			<< setw( 5) << EdgeList[e].node2
			<< setw( 6);
	
		switch (EdgeList[e].type)
		{
			case Y:		
				out << "Y"	<< "          "
					 << setw(20) << EdgeList[e].value;	
				break;
			case Z:		
				out << "Z"	<< "          "
				  	 << setw(20) << EdgeList[e].value;	
				break;
			case CC:	
				out << "CC";	
				out << "  ->" << setw(6)  
				     << EdgeList[EdgeList[e].pnum].name;
				break;
			case VC:	
				out << "VC";	
				out << "  ->" << setw(6)  
				     << EdgeList[EdgeList[e].pnum].name;
				break;
			case CS:	
				out << "CS";	
				out << "  <-" << setw(6)  
				     << EdgeList[EdgeList[e].pnum].name
				  	 << setw(20) << EdgeList[e].value;	
				break;
			case VS:	
				out << "VS";	
				out << "  <-" << setw(6)  
				     << EdgeList[EdgeList[e].pnum].name
				  	 << setw(20) << EdgeList[e].value;	
				break;
			default:
				break;
		}
		out << endl;
	}
}


/**	
	PrintPartialEdgeList()  prints an intermediate reduced graph edges.
	
	-- (newly added by G. Shi, 09-09-2010)
*/
void GPDD::PrintPartialEdgeList(EdgeInfo *edgelist, int length)
{
	cout << LINE_INFO << endl;
	
	int	num;

	if (!edgelist) {
		cout << LINE_INFO << "\n\tEmpty edge list." << endl;
	}

	for (int i = 0; i < length; i ++)
	{
		num = edgelist[i].edgenum;
		printf("\n\t[%2d] \"%6s\"\t (#%2d) : (%2d) ---- (%2d)  : paired with - %s", 
				i, 
				EdgeList[num].name,  	
				edgelist[i].edgenum,
				edgelist[i].node1,		// the current node number
				edgelist[i].node2, 
				(EdgeList[num].pname? EdgeList[num].pname : "(null)"));
	}
	cout << endl;
}


void GPDD::PrintGPDDStatistics()
{
	CountGPDDTerms();

	cout << LINE_INFO << "\n\tStatistics of GPDD:" << endl;
	
	cout << "\t--------------------------------------------------------" << endl;
//	cout << "\tNum of Minus Terms : " << GetCountOfMinusTerms() << endl;
//	cout << "\tNum of Plus Terms  : " << GetCountOfPlusTerms() << endl;
	cout << "\tNumber of GPDD nodes: " << GetTotalGPDDnodes() << endl;
	cout << "\tTerm Numbers  : " << scientific << GetNumTerms() << endl;
	cout << "\tCount of GPDD node creations: " << GetNodeCreationCount() << endl;
	cout << "\tNumber of Subgraphs : " << GetSubgraphCount() << endl;
	cout << "\tSubgraph Sharing: " << GetSubgraphSharingCount() << endl;	
	cout << "\tGPDD node Sharing: " << GetGPDDnodeSharingCount() << endl;
	cout << "\tNumber of GPDD nodes removed: " << GetNumNodesRemoved() << endl;
	cout << "\tNumber of GPDD nodes remaining: " << GetNumNodesRemaining() << endl;

#ifdef CALCULATE_S_COEFFICIENTS // statistics for s-expanding results
	//cout << "\tShort Power: " << pExpTpddRoot->pShort->sPowerList->size() << endl;
	//cout << "\tOpen Power: "  << pExpTpddRoot->pOpen->sPowerList->size() << endl;
	cout << "\tNumber of s-coefficient nodes: " << Get_sCoeffNodeCount() << endl;
	cout << "\tNumber of expanded GPDD nodes : " << GetTotalExpandedGPDDnodes() << endl;
#endif   // CALCULATE_S_COEFFICIENTS

}



void GPDD::PrintSymbols()
{
	cout << LINE_INFO << endl;
	
	if (NumSymbols == 0)
	{
		cout << "\nEmpty symbol list.";
		return;
	}
	
	cout << LINE_INFO << "\n\tNumSymbols = " << NumSymbols << endl;

	cout << setiosflags(ios::right)
		 << setw(10) << "Index"
		 << setw(10) << "Symbol"
		 << setw(10) << "IsSrc"
		 << setw(15) << "Value"  
		 << setw(10) << "Type"
		 << setw(15) << "Composite"
		 << endl;

	SymbNode	*pNode, *pNode2;
	int		n = 1;		// Count the index.

	for (n = 0; n < NumSymbols; n++)
	{
		pNode = LumpedSymbList[n];

		cout << setiosflags(ios::right)
			 << setw(10) << n
			 << setw(10) << pNode->symbname 
			 << setw(10) << (pNode->isSrc ? "Y" : "N")
			 << setw(15) << pNode->value
			 << setw(10) 
			 << EdgeTypeName[pNode->edge->type] << " : " 
			 << EdgeTypeName[pNode->isSrc ? pNode->pairedge->type : pNode->edge->type];
			 /*<< setw(10) << pNode->edge->node1
			 << setw(10) << pNode->edge->node2;*/

		if ((pNode2 = pNode->nextP))
		{
			cout << setw(15) << "(C)";
			cout << endl << "\t\t   +";
			// Print the lumped symbols
			while (pNode2)
			{
				cout << "->" << pNode2->symbname;
				pNode2 = pNode2->nextP;
			}
		}
		cout << endl;
	}
}
	
/**
  Binary Decision - Include Symbol
*/	
GPDDnode* GPDD::Reduction_INCLUDE_R(GPDDnode * pParent)
{
	GPDDnode * pNewNode, *pOpen, *pShort, *pShareNode;
	int isHash;

	//printf("\n[%s, %d] [Paranet] node(s%d)\n", __FILE__, __LINE__, pParent->GetSymbolIndex());

	pNewNode = CreateGPDDnode(pParent, SYMBOL_IN);

	//printf("\n[%s, %d] New [In] node(s%d)\n", __FILE__, __LINE__, pNewNode->GetSymbolIndex());
	NodeCreationCount ++;
	if (NodeCreationCount % 50000 == 0)
	{
		cout << LINE_INFO << "NodeCreationCount = " << NodeCreationCount << endl; 
	}

	if(pNewNode == VertexZero || pNewNode == VertexOne || !pNewNode)
		return pNewNode;
	
	pShareNode = ShareGPDDnode(pNewNode, isHash);
	if(isHash)
	{
	    return pShareNode;
	}

	//New node
	TotalGPDDnodes++;
	//cout<<LINE_INFO<<"TotalGPDDnodes = "<<TotalGPDDnodes<<endl;
	
	pShort = Reduction_INCLUDE_R(pNewNode);

	pOpen = Reduction_EXCLUDE_R(pNewNode);
	
	pNewNode->SetVertex_IN(pShort);
	pNewNode->SetVertex_OUT(pOpen);

	return pNewNode;
}

/**
  Binary Decision - Exclude Symbol
*/	
GPDDnode* GPDD::Reduction_EXCLUDE_R(GPDDnode * pParent)
{
	GPDDnode * pNewNode, *pOpen, *pShort, *pShareNode;
	int isHash;

	pNewNode = CreateGPDDnode(pParent, SYMBOL_OUT);

	//printf("\n[%s, %d] [Paranet] node(s%d)\n", __FILE__, __LINE__, pParent->GetSymbolIndex());
	//printf("\n[%s, %d] New [Out] node(s%d)\n", __FILE__, __LINE__, pNewNode->GetSymbolIndex());
	
	NodeCreationCount ++;
	if (NodeCreationCount % 50000 == 0)
	{
		cout << LINE_INFO << "NodeCreationCount = " << NodeCreationCount << endl; 
	}

	if (pNewNode == VertexZero || pNewNode == VertexOne || !pNewNode)
		return pNewNode;

	pShareNode = ShareGPDDnode(pNewNode, isHash);
	if (isHash)
	{
		return pShareNode;
	}
	
	//New node
	TotalGPDDnodes++;
	
	//cout<<LINE_INFO<<"TotalGPDDnodes = "<<TotalGPDDnodes<<endl;

	pShort = Reduction_INCLUDE_R(pNewNode);

	pOpen = Reduction_EXCLUDE_R(pNewNode);

	pNewNode->SetVertex_IN(pShort);
	pNewNode->SetVertex_OUT(pOpen);
	
	return pNewNode;
}

void GPDD::SetAllNodeLabels(HashTree ***ReduceHTree)
{
	cout << LINE_INFO << endl;

	int i;
	GPDDnode * p;

	//smaller the symb index is, higher the priority (order)
	//So here begin from the symbol with largest symbol index
	for (i = NumSymbols + 1; i > 1; i --)
	{
		p = ReduceList[i];	// This is an array of sublists (2D).
		while (p)
		{
			//SetNodeLabelR(p, ReduceHTree[i - 2], i - 2);
			SetNodeLabel(p, ReduceHTree[i - 2]);
			p = p->GetNextNode();		// scan the sublist
		}
	}	
}


//void GPDD::SetNodeLabelR(GPDDnode *p, HashTree **ReduceHTree, int level)
void GPDD::SetNodeLabel(GPDDnode *p, HashTree **ReduceHTree)
{
	int tagl, tagr, label_top;
	
	tagl = p->GetSign_IN() * p->GetVertex_IN()->GetTag()->label_top;		// multiplied by sign
	tagr = p->GetSign_OUT() * p->GetVertex_OUT()->GetTag()->label_top;

	/* cout << LINE_INFO << "\n\tSymb [#" << p->GetSymbolIndex() 
		<< "] with label triple = (" << p->GetTag()->label_top 
		<< ", " << tagl << ", " << tagr << ") " << endl; */

	//lable of VertexZero is 1
	if (abs(tagl) == 1 && abs(tagr) == 1)    // pShort == pOpen == VertexZero
	{
		p->GetTag()->SetTriple(1, 0, 0);	// equiv. to Zero-Suppression
		NumNodesRemoved ++;
		return;
	}

	GPDDnode_tag *pt = p->GetTag();
	GPDDnode_tag *qt;

	pt->SetTriple(0, tagl, tagr);

	int 	level = pt->HashId();

	//cout << LINE_INFO << "Get Hash ID = " << level << endl;
	
	//No entry of Symb with index level in the Reduce hash table yet
	if (!ReduceHTree[level])
	{
		ReduceHTree[level] = new HashTree();
		//cout << LINE_INFO << "Add new hash tree" << endl;
	}

	if ((qt = (GPDDnode_tag*) ReduceHTree[level]->Get(pt)))
	{
		label_top = qt->label_top;
		pt->SetTriple(label_top, tagl, tagr);
		NumNodesRemoved ++;
		//cout << LINE_INFO << "Find the existing node" << endl;
	}
	else
	{
		ReduceHTree[level]->Put(pt, *pt);
		
		//The unhashed node is saved to the NodeList vector.
		NodeList.push_back(p);
		
		NumNodesRemaining = NodeList.size();  // get the vector location for the label_top
		
		pt->SetTriple(NumNodesRemaining, tagl, tagr);
		
		p->Reference();
		//cout << LINE_INFO << "Added a new node to the hashtree" << endl;
	}

	/*cout << LINE_INFO << "\n\tSymb [#" << p->GetSymbolIndex() << "] is set label top = " 
		<< p->GetTag()->label_top << endl; */
	
	if (NumNodesRemaining % 10000 == 0)
	{
		cout << LINE_INFO << "Reconstructing the reduced GPDD ... " << NumNodesRemaining << endl;
	}
}

sCoeffNode *GPDD::Share_sCoeffNode(sCoeffNode *&p, int& tag)
{
#ifdef HASHSCOEFFNODE 
	int hashindex, level;
	sCoeffNode *q;
	level = p->pExpNode->SymbolIndex % TAG_PRIME1;
	hashindex = p->HashTreeID();
	 

	if(!uniqCoeffGPDD[level][hashindex])
	{
		uniqCoeffGPDD[level][hashindex] = new HashTree();
	}
	if((q = (sCoeffNode *) uniqCoeffGPDD[level][hashindex]->Get(p)))
	{
		delete p;
		p = 0;
	//	GPDDnodeSharingCount ++;
		tag = 1;
		return q;
	}
	else
	{
		uniqCoeffGPDD[level][hashindex]->Put(*p, p);
	}
	p->label = sCoeffNodeCount;

//	cout << sCoeffNodeCount << endl;
	sCoeffNodeCount ++;

	tag = 0;
	return p;	
#else
	sCoeffNodeCount ++;
	return p;
#endif

}

/**
	Function used by ExpandLumpedGPDD().
*/
//ExpGPDDnode *GPDD::Share_sExapndGPDDnodeS(ExpGPDDnode *p, int& tag)
// old name by W. Chen
ExpGPDDnode *GPDD::Share_ExpandedGPDDnode(ExpGPDDnode *p, int& tag)
{
	int hashindex, level;
	ExpGPDDnode *q;
	level = p->SymbolIndex % TAG_PRIME1;
	hashindex = p->HashTreeID();
	 
	if(!uniqExpGPDD[level][hashindex])
	{
		uniqExpGPDD[level][hashindex] = new HashTree();
	}
	if((q = (ExpGPDDnode *) uniqExpGPDD[level][hashindex]->Get(p)))
	{
		delete p;
		p = 0;
	//	GPDDnodeSharingCount ++;
		tag = 1;
		return q;
	}
	else
	{
		uniqExpGPDD[level][hashindex]->Put(p, *p);
	}
//	p->label = sCoeffNodeCount ++;

	TotalExpandedGPDDnodes ++;

//	cout << TotalExpandedGPDDnodes << endl;
	tag = 0;/**/
	return p;
}

/**
	If hashed, <hashed> returns 1; otherwise, it returns 0.
*/
GPDDnode *GPDD::ShareGPDDnode(GPDDnode *p, int &hashed)
{
	int 	hashindex, level;
	GPDDnode *q;
	level = p->GetSymbolIndex() % GPDD_PRIME_LEVEL;
	hashindex = p->HashTreeID();
	 
	if (!uniqGPDD[level][hashindex])
	{
		uniqGPDD[level][hashindex] = new HashTree();
	}
	
	if ((q = (GPDDnode *) uniqGPDD[level][hashindex]->Get(p)))
	{
		delete p;
		p = 0;
		GPDDnodeSharingCount ++;
		
		hashed = 1;
		return q;
	}
	else
	{
		uniqGPDD[level][hashindex]->Put(p, *p);
		
		hashed = 0;
		return p;
	}
}


Subgraph* GPDD::ShareSubgraph(Subgraph *p)
{
	int hashindex;
	Subgraph *q;


	hashindex = p->HashTreeID() % SUBGRAPH_PRIME;

	if(!uniqSubgraph[hashindex])
	{
		uniqSubgraph[hashindex] = new HashTree();
	}

	if((q = (Subgraph*) uniqSubgraph[hashindex]->Get(p)))
	{
	    
	    //printf("\n[%s : %d ] Sharing Sub-tree Pointer = 0x%x ... \n", __FILE__, __LINE__, q);

		delete p;
		p = 0;
		SubgraphSharingCount ++;
		return q;
	}
	else
	{
		uniqSubgraph[hashindex]->Put(p, *p);
	}

	return p;

}

/**
	The GPDD lumped nodes must be expanded (for the lumped symbols) by calling
	ExpandLumpedGPDD() before calling sExpandGPDD().
*/
//void GPDD::scoeffNew()		// Old name by W. Chen
void GPDD::sExpandGPDD()		// New name by G. Shi
{
	cout << LINE_INFO << endl;

	if (sExpandedFlag) {
		cout << LINE_INFO << "\n\tThe GPDD has been S-expanded." << endl;
		return;
	}
	
	if (!LumpedExpandedFlag) {	// The GPDD has not been lump-expanded yet !
		ExpandLumpedGPDD();
	}

	stack<ExpGPDDnode*> S;
	
	ExpGPDDnode *pCurNode;

	sCoeffNodeCount = 2;

	// sCoeffOne and sCoeffZero will be deleted
	// by the destructors of ExpOne and ExpZero.
	sCoeffOne = new sCoeffNode;
	sCoeffZero = new sCoeffNode;

	sCoeffOne->pExpNode = ExpOne;
	sCoeffOne->Vertex_IN = NULL;
	sCoeffOne->Vertex_OUT = NULL;
	sCoeffOne->label = 1; 
	sCoeffOne->Value = 1.0;

	sCoeffZero->pExpNode = ExpZero;
	sCoeffZero->Vertex_IN = NULL;
	sCoeffZero->Vertex_OUT = NULL;
	sCoeffZero->label = 0;
	sCoeffZero->Value = 0.0;


	ExpOne->sCoeffList = new vector<sCoeffNode*>;
	ExpOne->sPowerList = new vector<int>;
	ExpOne->sCoeffList->push_back(sCoeffOne);
	ExpOne->sPowerList->push_back(0);

	ExpZero->sCoeffList = new vector<sCoeffNode*>;
	ExpZero->sPowerList = new vector<int>;
	
	ExpZero->sCoeffList->push_back(sCoeffZero);	// Keep it empty
	ExpZero->sPowerList->push_back(0);			// Keep it empty  (G. Shi)
	// Because it just adds 0 to the ancester nodes.

	ExpOne->Mark();
	ExpZero->Mark();

	S.push(ExpGPDDroot);

	while(!S.empty())
	{
		pCurNode = S.top();

		if(!pCurNode->Vertex_IN->sCoeffList)
		{
			S.push(pCurNode->Vertex_IN);
		}
		else if(!pCurNode->Vertex_OUT->sCoeffList)
		{
			S.push(pCurNode->Vertex_OUT);
		}
		else
		{
			CreateCoeffList(pCurNode);	// Expand one GPDDnode
			S.pop();
		}
	}

	sExpandedFlag = true;	// have to expand the s-polynomial only once
	
}

/**
	This function is used in the GPDD deconstructor for cleaning up
	the s-expanded GPDD nodes.

	Since the lump-expanded GPDD nodes are not in the unqiue table, 
	we should release them manually.
*/
void GPDD::CleanExpandedGPDDnodes()
{	
	ExpOne->visit_flag	= true;
	ExpZero->visit_flag = true;
	
	CleanExpandedGPDDnodes_R(ExpGPDDroot);
}

void GPDD::CleanExpandedGPDDnodes_R(ExpGPDDnode *&p)
{
	if (p->visit_flag)
		return;
	
	char eType = LumpedSymbList[p->SymbolIndex]->symbname[0];
	ExpGPDDnode 	*q, *pcur;

	pcur = p;

	if (eType == '^')
	{
		q = p->Vertex_OUT;

		while (q->pGPDDnode == NULL && q->Vertex_OUT)
		{
			p = q;
			q = q->Vertex_OUT;

			delete p;
			p = NULL;
		}
		
		pcur->Vertex_OUT = q;
	}
	pcur->visit_flag = true;

	
	CleanExpandedGPDDnodes_R(pcur->Vertex_IN);
	CleanExpandedGPDDnodes_R(pcur->Vertex_OUT);
}

/**
	An s-expanded polynomial is something like this :
		P_a(s) = a[n] s^n + a[n-1] s^{n-1} + ... + a[1] s + a[0]
	Each coefficient a[k] is a BDD, which is a representation of SOP.
	By running this function, the number of terms in each coefficient a[k]
	is calculated by traversing the corresponding BDD.

	(G. Shi, 09-08-2010)
*/
void GPDD::CountNumTermsForAllCoeff()
{
	cout << LINE_INFO << endl;
	
	sCoeffOne->NumTerms = 1;
	sCoeffZero->NumTerms = 0;

	for (unsigned i = 0; i < ExpGPDDroot->Vertex_IN->sCoeffList->size(); i ++)
	{
		CountNumTermsOneCoeff_R((*(ExpGPDDroot->Vertex_IN->sCoeffList))[i]);
	}

	for (unsigned i = 0; i < ExpGPDDroot->Vertex_OUT->sCoeffList->size(); i ++)
	{
		CountNumTermsOneCoeff_R((*(ExpGPDDroot->Vertex_OUT->sCoeffList))[i]);
	}
}

void GPDD::CountNumTermsOneCoeff_R(sCoeffNode* p)
{
	if (p->NumTerms != -1)
		return;

	CountNumTermsOneCoeff_R(p->Vertex_IN);
	CountNumTermsOneCoeff_R(p->Vertex_OUT);

	p->NumTerms = p->Vertex_IN->NumTerms + p->Vertex_OUT->NumTerms;
}


void GPDD::RunSymbolStatistics()
{
	cout << LINE_INFO << endl;
	cout << "\n\tStatistics of GPDD symbols:" << endl;
	
	stack<ExpGPDDnode*> *S = new stack<ExpGPDDnode*>;
	
///	UnMarkAllNodes(ExpGPDDroot, 2);

	TraversalMark = ExpGPDDroot->GetMark();
	FlipTraversalMark();

	ExpGPDDnode* pTopNode;

	ExpOne->SetMark(TraversalMark);
	ExpZero->SetMark(TraversalMark);

	ExpOne->NumTerms = 1;
	ExpZero->NumTerms = 0;

	ExpOne->SymbCountArray = new double[TotalNumSymbols];
	ExpZero->SymbCountArray = new double[TotalNumSymbols];

	memset(ExpOne->SymbCountArray, 0, sizeof(double) * TotalNumSymbols);
	memset(ExpZero->SymbCountArray, 0, sizeof(double) * TotalNumSymbols);

	S->push(ExpGPDDroot);

	while (!S->empty())
	{
		pTopNode = S->top();

		if (pTopNode->GetMark() == TraversalMark)
		{
			S->pop();
		}
		else
		{
			if (!pTopNode->GetVertex_IN()->GetMark() == TraversalMark)
			{
				S->push(pTopNode->GetVertex_IN());
			}
			else if (!pTopNode->GetVertex_OUT()->GetMark() == TraversalMark)
			{
				S->push(pTopNode->GetVertex_OUT());
			}
			else
			{
				pTopNode->SetMark(TraversalMark);
				
				pTopNode->SymbCountArray = new double[TotalNumSymbols];
				memset(pTopNode->SymbCountArray, 0, sizeof(double) * TotalNumSymbols);

				pTopNode->NumTerms = pTopNode->GetVertex_IN()->NumTerms + 
										pTopNode->GetVertex_OUT()->NumTerms;

				for (int i = 0; i < TotalNumSymbols; i ++)
				{
					pTopNode->SymbCountArray[i] = pTopNode->GetVertex_IN()->SymbCountArray[i] 
						+ pTopNode->GetVertex_OUT()->SymbCountArray[i]; 
				}
				pTopNode->SymbCountArray[pTopNode->SymbolIndex] = 
							pTopNode->GetVertex_IN()->NumTerms;
			}
		}
	}

	cout << "\n    " << setw(20) << "Symbol Name" << setw(20) << "Symbol Count" << endl;
	
	for (int i = 0; i < TotalNumSymbols; i ++)
	{
		cout << "    " << setw(20) << FullSymbList[i]->symbname  
			 << setw(20) << ExpGPDDroot->SymbCountArray[i] << endl;
	}
}


void GPDD::MarkAllNodes(BDDnode* p, int type)
{
	stack<BDDnode*> *S = new stack<BDDnode*>;

	GPDDnode* pTopNode;
	ExpGPDDnode* peTopNode;

	S->push(p);
	
	if (type == 1)
	{
		VertexOne->SetMark(TraversalMark);
		VertexZero->SetMark(TraversalMark);
		
		while (!S->empty())
		{
			pTopNode = (GPDDnode*) S->top();
			if (pTopNode->GetMark() == TraversalMark)
			{
				S->pop();
			}
			else
			{
				if (pTopNode->GetVertex_IN()->GetMark() != TraversalMark)
					S->push(pTopNode->GetVertex_IN());
				else if (pTopNode->GetVertex_OUT()->GetMark() != TraversalMark)
					S->push(pTopNode->GetVertex_OUT());
				else
				{
					pTopNode->SetMark(TraversalMark);
				}
			}
		}
	}
	else
	{	
		ExpOne->SetMark(TraversalMark);
		ExpZero->SetMark(TraversalMark);
		
		while (!S->empty())
		{
			peTopNode = (ExpGPDDnode*)S->top();
			if (peTopNode->GetMark() == TraversalMark)
			{
				S->pop();
			}
			else
			{
				if (peTopNode->GetVertex_IN()->GetMark() != TraversalMark)
					S->push(peTopNode->GetVertex_IN());
				else if (peTopNode->GetVertex_OUT()->GetMark() != TraversalMark)
					S->push(peTopNode->GetVertex_OUT());
				else
				{
					peTopNode->SetMark(TraversalMark);
				}
			}
		}	
	}
}

/**
	** SENSITIVITY **
	Look up the index of a symbol given its name, meanwhile set the sensitivity index.

	Returns:  (index) if found;
	             (NOPE) if failed.
*/
int GPDD::LookupSymbolIndex(const char * name)
{
	int i = 0;
	
	if (!pGPDDroot || pGPDDroot->GetSymbolIndex() < 0)
	{
		cout << LINE_INFO << "there is no root.." << endl;
		return (NOPE);
	}

	SymbNode *pParaNode;
	
	for (i = 0; i < NumSymbols; i++)
	{
//		cout << i << "->" << LumpedSymbList[i]->symbname <<endl;
		
		if (LumpedSymbList[i]->symbname[0] !='^')
		{
			if (strcmp(LumpedSymbList[i]->symbname, name) == 0)
			{
				SetSensitivityIndex(i);		// Can set it after looking the index up (G.Shi)

				cout << LINE_INFO << " Found the symbol \"" << name << "\"" << endl;
				return (i);
			}
		}
		else
		{
			pParaNode = LumpedSymbList[i]->nextP;
			while (pParaNode)
			{
				if (strcmp(pParaNode->symbname, name) == 0)
				{
					SetSensitivityIndex(i);
					
					cout << LINE_INFO << " Found the symbol \"" << name << "\"" << endl;
					return (i);
				}
				pParaNode = pParaNode->nextP;
			}
		}
	
	}
//	cout <<  "the index is : " << i << endl;

	cout << LINE_INFO << "Did not find the symbol \"" << name << "\"" << endl;
	return (NOPE);
}



// Change the value of a specific component
// Li Ji
// 2008-10-29 22:46

/**
	Only sets a new value to the symbol, the old value is erased.
*/
//double GPDD::change_value(const string &symbname, double newvalue)
int  GPDD::SetNewValue(const string &symbname, double newvalue)
{
	
	if (NumSymbols == 0)
	{
		cout << LINE_INFO << "\nEmpty symbol list.\n" << endl;
		return (STOP); // an error value
	}
	
	for (int n=0; n < NumSymbols; ++n)
	{
		SymbNode *pNode = LumpedSymbList[n];

		for ( ; pNode; pNode = pNode->nextP)
		{
			string name(pNode->symbname);
			if (name == symbname)
			{
				cout << LINE_INFO << "\n\t"
					<< name << " : old value (" << pNode->value 
					<< ") --> new value ( " << newvalue << ") " << endl;
				pNode->value = newvalue;
				return (OK);
			}
		}
	}
	
	cout << LINE_INFO << "Symbol \" " << symbname << "\" not found." << endl;
	
	return (STOP);
}



/**
	Change the specified symbol value by percentage
*/
int  GPDD::VariationByPercentage(const string &symbname, double percent)
{
	
	if (NumSymbols == 0)
	{
		cout << LINE_INFO << "\nEmpty symbol list.\n" << endl;
		return (STOP); // an error value
	}
	
	for (int n=0; n < NumSymbols; ++n)
	{
		SymbNode *pNode = LumpedSymbList[n];

		for ( ; pNode; pNode = pNode->nextP)
		{
			string name(pNode->symbname);
			if (name == symbname)
			{
				cout << LINE_INFO << "\n\t"
					<< name << " : old value (" << pNode->value 
					<< ") --> percent of ( " << (percent*100) << "%) " << endl;
				
				pNode->value *= (1 + percent);
				return (OK);
			}
		}
	}
	
	cout << LINE_INFO << "Symbol \" " << symbname << "\" not found." << endl;
	
	return (STOP);
}

/**
	Set new value while get the old value.
	(Added by G. Shi)
*/
int  GPDD::SetNewValue(const string &symbname, double &oldvalue, double newvalue)
{
	
	if (NumSymbols == 0)
	{
		cout << LINE_INFO << "\nEmpty symbol list.\n" << endl;
		return (NOPE); // an error value
	}
	
	for (int n=0; n < NumSymbols; ++n)
	{
		SymbNode *pNode = LumpedSymbList[n];

		for ( ; pNode; pNode = pNode->nextP)
		{
			string name(pNode->symbname);
			if (name == symbname)
			{
				oldvalue = pNode->value;
				cout << LINE_INFO << "\n\t" 
					<< name << " : old value (" << pNode->value 
					<< ") --> new value ( " << newvalue << ") " << endl;
				pNode->value = newvalue;
				return (OK);
			}
		}
	}
	
	cout << LINE_INFO << "Symbol \" " << symbname << "\" not found." << endl;
	return (NOPE);
}


void GPDD::SetPairing()
{
	int e1;

	for (int e = 0; e < NumEdges; e++)
	{
		string	pname;
		if ( EdgeList[e].type > 1)		// For those dependent types.
		{
			pname = EdgeList[e].pname;
			e1 = 0;
			// Look for the pairing edge
			while ( pname.compare(EdgeList[e1++].name) );
			// Make a link.
			EdgeList[e].pnum = --e1;
		}
	}
}

void GPDD::SortSymbList()
{
#if 0	
	int i, n1, n2, j, tmp;

	int* sdegree = new int[N];
	int* rdegree = new int[N];
	memset(sdegree, 0, N * sizeof(int));
	memset(sdegree, 0, N * sizeof(int));
	stag *s_array = new stag[NumSymbols];
	qsort(LumpedSymbList + 1, NumSymbols - 1, sizeof(SymbNode*), compare);

	int *marked = new int[2 * NumSymbols];
	
	SymbNode **stemp = new SymbNode*[NumSymbols];

	memset(marked, 0, 2 * NumSymbols * sizeof(int));

	delete[] rdegree;
	rdegree = 0;
	delete[] sdegree;
	sdegree = 0;
	delete[] s_array;
	s_array = 0;
	delete[] stemp;
	stemp = 0;
	delete[] marked;/**/
	marked = 0;/**/
#endif	
}


void GPDD::ReferenceAllNodes()
{
	ResetAllReferences(pGPDDroot, 1);

//	UnMarkAllNodes(pGPDDroot, 1);

	TraversalMark = pGPDDroot->GetMark();
	FlipTraversalMark();

	if (!pGPDDroot)
		return;
	
	Reference_R(pGPDDroot);

	pGPDDroot->SetMark(TraversalMark);
}


void GPDD::Reference_R(GPDDnode *p)
{
	if (p->GetMark() == TraversalMark)
		return;
	
	p->SetMark(TraversalMark);

	GPDDnode *pNode;
	
	pNode = p->GetVertex_IN();

	if(pNode != VertexOne && pNode != VertexZero)
	{
		pNode->Reference();
		Reference_R(pNode);
	}

	pNode = p->GetVertex_OUT();
	if(pNode != VertexOne && pNode != VertexZero)
	{
		pNode->Reference();
		Reference_R(pNode);
	}

}

void GPDD::ResetAllReferences(BDDnode* p, int type)
{

///	UnMarkAllNodes(pGPDDroot, 1);

	stack<BDDnode*> *S = new stack<BDDnode*>;

	GPDDnode* pTopNode;
	ExpGPDDnode* peTopNode;

	TraversalMark = p->GetMark();
	FlipTraversalMark();

	S->push(p);

	if (type == 1)
	{
		VertexOne->SetMark(TraversalMark);
		VertexZero->SetMark(TraversalMark);
		
		while (!S->empty())
		{
			pTopNode = (GPDDnode*)S->top();
			
			if (pTopNode->GetMark() == TraversalMark)
			{
				S->pop();
			}
			else
			{
				if (!pTopNode->GetVertex_IN()->GetMark() == TraversalMark)
					S->push(pTopNode->GetVertex_IN());
				else if (!pTopNode->GetVertex_OUT()->GetMark() == TraversalMark)
					S->push(pTopNode->GetVertex_OUT());
				else
				{
					pTopNode->ResetReference();
					pTopNode->SetMark(TraversalMark);
				}
			}
		}
	}
	else
	{
		ExpOne->SetMark(TraversalMark);
		ExpZero->SetMark(TraversalMark);
		
		while (!S->empty())
		{
			peTopNode = (ExpGPDDnode*)S->top();
			
			if (peTopNode->GetMark() == TraversalMark)
			{
				S->pop();
			}
			else
			{
				if (!peTopNode->GetVertex_IN()->GetMark() == TraversalMark)
					S->push(peTopNode->GetVertex_IN());
				else if (!peTopNode->GetVertex_OUT()->GetMark() == TraversalMark)
					S->push(peTopNode->GetVertex_OUT());
				else
				{
					peTopNode->ResetReference();
					peTopNode->SetMark(TraversalMark);
				}
			}
		}	
	}
}


void GPDD::ResetReference_R(GPDDnode *p)
{
	if (p->GetMark() == TraversalMark || p == VertexOne || p == VertexZero)
		return;
	
	p->SetMark(TraversalMark);

	p->ResetReference();

	ResetReference_R(p->GetVertex_IN());
	ResetReference_R(p->GetVertex_OUT());
}


/**
	This function will be deprecated. (G. Shi, 09-16-2010)
	
	<type = 1> works on GPDDnode;
	<type = 0> works on ExpGPDDnode.
*/	
void GPDD::UnMarkAllNodes(BDDnode* p, int type)
{
//	cout << LINE_INFO << endl;

	stack<BDDnode*> *S = new stack<BDDnode*>;

	GPDDnode* pTopNode;
	ExpGPDDnode* peTopNode;

	/* 
		Mark the root ! Otherwise, unmarking would not procede.
		It is important for repeated GPDD traversals. (BUG! G. Shi, 09-06-2010)
	*/
	p->Mark();		/********/
	
	S->push(p);

	if (type == 1)		// Work on GPDDnodes
	{
		VertexOne->UnMark();
		VertexZero->UnMark();
		while (!S->empty())
		{
			pTopNode = (GPDDnode*)S->top();
			//printf("pTopNode = 0x%x\n", pTopNode);

			if (!pTopNode->GetMark())
			{
				S->pop();
			}
			else
			{
				if (pTopNode->GetVertex_IN()->GetMark())
					S->push(pTopNode->GetVertex_IN());
				else if (pTopNode->GetVertex_OUT()->GetMark())
					S->push(pTopNode->GetVertex_OUT());
				else
				{
					pTopNode->UnMark();
				}
			}
		}
	}	
	else		// Work on ExpGPDDnodes
	{
		ExpOne->UnMark();
		ExpZero->UnMark();
		while (!S->empty())
		{
			peTopNode = (ExpGPDDnode*)S->top();
			if (!peTopNode->GetMark())
			{
				S->pop();
			}
			else
			{
				if (peTopNode->GetVertex_IN()->GetMark())
					S->push(peTopNode->GetVertex_IN());
				else if (peTopNode->GetVertex_OUT()->GetMark())
					S->push(peTopNode->GetVertex_OUT());
				else
				{
					peTopNode->UnMark();
				}
			}
		}
	}
}

	
void GPDD::UnMarkNode_R(GPDDnode *p)
{
	if (!p->GetMark())
		return;

	p->UnMark();

	if(p->GetVertex_IN())
	{
		UnMarkNode_R(p->GetVertex_IN());
	}
	if(p->GetVertex_OUT())
	{
		UnMarkNode_R(p->GetVertex_OUT());
	}
}


