Add files via upload

1 files changed, 1002 insertions, 0 deletions

diff --git a/sourcecodes/k-best/src/Model.h b/sourcecodes/k-best/src/Model.h
new file mode 100644
index 00000000..1ba7a0dd
--- /dev/null
+++ b/sourcecodes/k-best/src/Model.h
@@ -0,0 +1,1002 @@
+#ifndef MODEL_H
+#define MODEL_H
+
+#include<stdio.h>
+#include<stdlib.h>
+
+#include<iostream>
+#include<vector>
+#include<deque>
+#include<string>
+#include<fstream>
+#include<math.h>
+
+#include"Arguments.h"
+
+using namespace std;
+
+#include"UpdateHR2.h"
+
+#define MAX_COUNT 1000
+#define MAX_RQ 8*10000
+static double lr[MAX_RQ][MAX_COUNT]; 
+static double lr2[MAX_RQ][MAX_COUNT]; 
+
+//ofstream totalOf1("debugModelTotal1.txt");
+
+double log_gammaratio(int a, int z){
+	if(z >= MAX_COUNT  || a >= MAX_RQ){
+//		if(z >= MAX_COUNT){
+//			cerr << " Counts larger than " << MAX_COUNT << " occured. " << endl;
+//			cerr << " Increase MAX_COUNT. Exit now.\n"; 
+//			
+//		}
+//		else if(a >= MAX_RQ ){
+//			cerr << " a larger than " << MAX_RQ << " occured. " << endl;
+//			cerr << " Increase MAX_RQ. Exit now.\n";  
+//		} 
+//		
+//		exit(1); 		
+
+		double log_gam_res = 0;
+		//seperate 2 cases:
+		//if only z >= MAX_COUNT but a < MAX_RQ
+		if(a < MAX_RQ){
+			//can compute it faster based on lr
+			log_gam_res += lr[a][MAX_COUNT - 1];
+			for(int k = MAX_COUNT - 1; k <= z-1; k++){
+				log_gam_res += log(1.0/a + k);
+			}
+			
+		}
+		//o.w.: start from the beginning
+		else{
+			for(int k = 0; k <= z-1; k++){
+				log_gam_res += log(1.0/a + k);
+			}
+		}
+
+//		totalOf1 << "log_gam_res = "  << log_gam_res << endl << endl;
+		return log_gam_res;
+
+	}
+	else if (z < MAX_COUNT && a > 0){
+		return lr[a][z];
+	}
+	else if (z < MAX_COUNT && a < 0){
+		return lr2[-a][z];
+	}
+	else{
+		cerr << "Error: a = 0" << endl;
+		exit(1); 
+	}
+}
+
+
+double log_gammaratio(double z){
+	return 0.50 * log(M_PI) + (z - 0.50) * log(z) - z + 1.0 / (12 * z);  
+}
+
+template<class T> void printvec(ostream& f, const vector<T>& v){
+	int size= v.size();
+	int i;
+	for(i=0; i<size; i++){
+		f<<v.at(i)<<" ";
+	}
+	f<<endl;
+}
+template<class T> void printvecs(ostream& f, const vector<vector<T> >& v){
+	int size= v.size();
+	int i;
+	for(i=0; i<size; i++){
+		printvec(f, v.at(i));
+	}
+}
+
+
+void print_nodes(ostream& f, int *V, int n){
+	for (int j = 0; j < n; j ++){
+		f<<" "<<V[j];
+	}
+}
+
+void print_nodes(ostream& f, int T, int *V, int n){
+	int j = 0; 
+	while (T > 0){           
+		if (T & 1){
+			f<<" "<<V[j];
+		}
+		j ++;
+		
+		T >>= 1; 
+	}
+}
+void print_set(ostream& f, int T){
+	int j = 0; 
+	while (T > 0){
+		if (T & 1){
+			f<<" "<<j;
+		}
+		j ++;
+		T >>= 1;
+	}
+}
+
+
+// Reads and stores data.
+class Data{
+public:
+	Data(){}
+	~Data(){}
+	void init(){
+		read_data();
+		downcode();
+	}
+	void read_data(){
+		ifstream ifs(Arguments::datafile, ios::in);
+		if (!ifs){
+			fprintf(stderr, " Cannot read file %s.\n", Arguments::datafile);
+			exit(1);
+		}
+		fprintf(stderr, " Reading file %s...\n", Arguments::datafile);
+		
+		char* buffer= new char[10000];
+		ifs.getline(buffer, 10000);
+		char* pch= strtok(buffer,", \t");
+		while (pch != NULL){
+			string tempstring(pch);
+			//cerr << "Node = "<< tempstring << endl;
+			heads.push_back(tempstring);
+			pch = strtok(NULL, ", \t");
+        }
+		numattributes = heads.size();
+	
+		fprintf(stderr, " Heading read: %d attributes.\n", numattributes);
+	
+		dm.clear();
+		vector<int> temp;
+		while (true){
+			temp.clear();
+			ifs.getline(buffer, 10000);
+			pch = strtok(buffer,", \t");
+			while (pch != NULL){
+				temp.push_back(atoi(pch));
+				pch = strtok(NULL, ", \t");
+			}
+			if ((int)temp.size() != numattributes) break;
+			
+			dm.push_back(temp);
+			
+			if ((int)dm.size() >= atoi(Arguments::maxnumrecords)) break;
+		}
+		numrecords = dm.size();
+		fprintf(stderr, " Data read: %d lines.\n", numrecords); 
+		delete [] buffer;
+	}
+	void downcode(){
+		//fprintf(stderr, "downcode() starts. \n"); 
+		
+		int inuse[4096][256];
+		for (int v = 0; v < 255; v ++){
+			for (int i = 0; i < numattributes; i ++){
+				inuse[v][i] = -1;
+			}
+		}
+		//Purpose: For each i attribute, 
+		for (int t = 0; t < numrecords; t ++){
+			for (int i = 0; i < numattributes; i ++){
+				inuse[dm[t][i]][i] = 1;
+			}
+		}
+
+		//ofstream of1("debugDowncode.txt");
+		//of1 << "\n inuse check" << endl;
+		for (int v = 0; v < 30; v ++){			
+
+			
+//			for (int i = 0; i < numattributes; i++){
+//				of1 << " inuse[" << v << "][" << i << "] =" << inuse[v][i] << endl;
+//			}
+//			of1 << "\n";
+			
+
+         	/*
+         	FILE * fp1 = fopen("debug1.txt", "w");
+          	fprintf(fp1, "\n inuse check\n");
+          	for (int i = 0; i < numattributes; i++){
+				fprintf(fp1, " inuse[%d][%d] = %d\n", v, i, inuse[v][i]);
+			}
+                  fclose(fp1);
+ 			*/
+                  
+
+			/*
+      		fprintf(stderr, "\n inuse check\n");
+           	for (int i = 0; i < numattributes; i++){    	
+				fprintf(stderr, " inuse[%d][%d]=%d\n", v, i, inuse[v][i]);
+          	}
+			fprintf(stderr, "\n");
+			*/
+		}
+
+		arities.clear();
+
+		maxarity = 0;
+
+		for (int i = 0; i < numattributes; i ++){
+			int count = 0;
+			for (int v = 0; v < 255; v ++){
+				if (inuse[v][i] == 1){
+					inuse[v][i] = count ++;
+				}
+			}
+			arities.push_back(count);	
+			if (count > maxarity) maxarity = count;
+		}
+
+//      of1 << "maxarity = " << maxarity << endl;
+//		printvec(of1, arities);
+		
+		for (int t = 0; t < numrecords; t ++){
+			for (int i = 0; i < numattributes; i ++){
+				dm[t][i] = inuse[dm[t][i]][i];
+
+			}
+		}
+	}//end downcode()
+	
+	
+	int get_index(string s){
+		int i = 0;
+		while (i < (int)heads.size() && heads[i] != s){ i ++; }
+		return i;
+	}
+	void print_data(ostream & f){
+		for (int i = 0; i < numattributes; i ++){
+			f << " " << heads[i];
+		}
+		f << endl; 
+		for (int t = 0; t < numrecords; t ++){
+			for (int i = 0; i < numattributes; i ++){
+				f << " " << dm[t][i];
+			}
+			f << endl; 
+		}
+	}
+		
+	vector<string> heads;
+
+	vector< vector<int> > dm;
+
+	vector<int> arities;
+
+	int numattributes;
+
+	int numrecords;
+
+	int maxarity;
+};
+
+
+// NOTE: We assume that children come first, that is, 
+// V[0] is the grandest child and V[n-1] is the grandest parent.
+// Layers respectively from 0 to numlayers.
+// Note that this is reverse to the input ordering.
+class Layering{
+public:
+	Layering(){}
+	~Layering(){ delete [] V;}
+
+	void init(Data & data){
+		set_layers(data);
+		
+		//print_layers();
+	
+	}
+	void set_layers(Data & data){
+		
+//		ofstream of1("debugSet_layers.txt");
+	
+		V = new int[data.numattributes];
+		int j = 0;
+		cnh.clear(); nh.clear(); cnh.clear();
+		
+		if (Arguments::layeringfile[0] == '%'){
+			numlayers = 1;
+			int layersize = data.numattributes;			
+			nh.push_back(layersize);
+			cnh.push_back(j);
+			
+			for (int l = 0; l < layersize; l ++){
+				int i = j;
+				V[j++] = i;
+			}
+			
+			edgeswithin.clear();
+			edgeswithin.push_back(true);
+        	
+			return;
+		}
+	
+		parselayeringfile(Arguments::layeringfile, true);
+		numlayers = nodesinlayers.size();
+	
+		for (int h = numlayers-1; h >= 0; h --){
+			int layersize = nodesinlayers[h].size();
+			nh.push_back(layersize);
+			cnh.push_back(j);
+			for (int l = 0; l < layersize; l ++){
+				int i = data.get_index(nodesinlayers[h][l]);
+				V[j++] = i;
+			}
+		}
+	
+	}
+	void print_layers(){
+		for (int h = 0; h < numlayers; h ++){
+			cerr<<"Layer "<<h<<":";
+			for (int j = 0; j < nh[h]; j ++){
+				cerr<<" "<<V[cnh[h] + j];
+			}
+			cerr<<endl;
+		}
+	}
+	bool edges_within(int h){
+		return edgeswithin[numlayers - h - 1];
+	}	
+	void parselayeringfile(char *layerfilename, bool verbose){
+	
+		ifstream ifs(layerfilename, ios::in);
+		if( !ifs){
+			cerr<<"ERROR: The file '"<<layerfilename<<"' could not be opened\n";
+			cerr<<"Exiting\n";
+			exit(1);
+		}
+		if(verbose){
+			cerr<<"Reading layer file: '"<<layerfilename<<"'\n";
+		}
+
+		char* buffer= new char[10000];
+		char* pch;
+		int i;
+		char tempchar;
+		while(true){
+			ifs.getline(buffer,10000);
+			pch = strtok(buffer, ", \t");
+			vector<string> tempstringvec;	//will hold the tokenized line
+			while(pch != NULL){
+				string tempstring(pch);
+				tempstringvec.push_back(tempstring);
+				pch = strtok(NULL, ", \t");
+			}
+
+			if(tempstringvec.size() < 3){
+				cerr<<"   ERROR: in reading file '"<<layerfilename<<"'"<<endl;
+				cerr<<"          each line must be of the form:"<<endl;
+				cerr<<"          layername\t1/0\tn1 ... nn"<<endl;
+				cerr<<"          EXITING"<<endl;
+				exit(1);
+			}
+
+			layernames.push_back(tempstringvec.at(0));
+			if(tempstringvec.at(1) == "0"){
+				edgeswithin.push_back(0);
+			}
+			else if( tempstringvec.at(1) == "1"){
+				edgeswithin.push_back(1);
+			}
+			else{
+				cerr<<"   ERROR: in reading file '"<<layerfilename<<"'"<<endl;
+				cerr<<"          the values in the second column can either be '1' or '0'"<<endl;
+				cerr<<"          EXITING"<<endl;
+				exit(1);
+			}
+
+			vector<string> tempnodesinlayer;
+			for(i = 2; i < (int)tempstringvec.size();i ++){
+				tempnodesinlayer.push_back(tempstringvec.at(i));
+			}
+			nodesinlayers.push_back(tempnodesinlayer);
+
+			ifs>>tempchar;
+			ifs.putback(tempchar);
+			if(ifs.eof()){
+				break;
+			}
+		}
+		delete [] buffer;
+		buffer = 0;
+		pch = 0;
+		cerr.width(13);
+		cerr<<"Layer Name";
+		cerr.width(24);
+		cerr<<"Edges Within Allowed";
+		cerr.width(14);
+		cerr<<"Node Names"<<endl;
+		for(i = 0; i < (int)layernames.size(); i ++){
+			cerr.width(13);
+			cerr<<layernames.at(i);
+			cerr.width(15);
+			if(edgeswithin.at(i)){
+				cerr<<"true";
+			}
+			else{
+				cerr<<"false";
+			}
+			cerr<<"             ";
+			printvec(cerr, nodesinlayers.at(i));
+		}
+	}
+
+
+	vector<string> layernames; 
+	vector<bool> edgeswithin;
+	vector<vector<string> > nodesinlayers;
+	int *V;
+	int numlayers;
+	vector<int> nh;
+	vector<int> cnh;
+};
+
+
+
+// Blocks for computing sufficient statistics.
+class Tnode{
+public:
+	double pseudocount;
+	short int arity;
+	short int count;
+	short int depth;
+	Tnode** children;
+	
+	Tnode(const int a, const int d){
+		pseudocount = 1;
+		arity = a;
+		count = 0;
+		depth = d;
+		children = new Tnode*[arity];
+		for(int v = 0; v < arity; v ++){
+			children[v] = NULL;
+		}
+	}
+	~Tnode(){
+		for(int v = 0; v < arity; v ++){
+			if (children[v] != NULL){
+				delete children[v];
+			}
+		}
+		delete [] children;
+		//free_memory();
+	}
+	void free_memory(){
+		//cerr<<" Deleting d="<<depth<<" c="<<count<<endl;
+		
+		for(int v = 0; v < arity; v ++){
+			if (children[v] != NULL){
+				children[v]->free_memory();
+				delete children[v];
+			}
+		}
+		//delete [] children;
+	}
+	// Note: Actually S is not needed.
+	void update(vector<int> & S, vector<int> & u){
+		if (u.size() > 0){	
+			//int j = S[u.size()-1];			
+			int v = u[u.size()-1]; 
+			
+			u.pop_back();
+			if (children[v] == NULL){
+				//cerr<<" Creating d="<<depth<<" c="<<count<<endl;
+				
+				children[v] = new Tnode(arity, depth+1);
+			}
+			children[v]->update(S, u);
+		}
+		count ++;
+	}
+	void print(ostream & f){
+		f << "(" << count << ")" << endl;
+		for (int v = 0; v < arity; v ++){
+			if (children[v] != NULL){
+				for (int d = 0; d < depth; d ++) f << " "; 
+				f << " " << v << " ";
+				children[v]->print(f);
+			}
+		}
+	}
+	
+	double evaluate(int d, int a, int r){
+		double sum = 0;
+		
+		//if no parent, directly log[gamma(1+count) / gamma(1)  ] (1/a = 1/1 = 1)
+		if (depth == d){
+				sum = log_gammaratio(a, count);
+				//cerr << "count1: " << count << endl;
+
+		}
+		else if (depth <= d - 1){
+			for (int v = 0; v < arity; v ++){
+				if (children[v] != NULL){
+					double value = children[v]->evaluate(d, a, r);
+					sum += value;
+				}
+			}
+			if (depth == d - 1){
+				sum -= log_gammaratio(r, count);
+				//cerr << "count2: " << count << endl;
+
+			}
+		}
+		return sum;
+	}
+	
+
+	double evaluateHR(int d, int alpha_v_pa, int r){
+		//cerr << "start: root.evaluateHR(d ="
+		//<< d << ", alpha_v_pa=" << alpha_v_pa << ", r=" << r << ")" << endl;
+		double sum = 0;
+		
+		//if no parent, directly log[gamma(1+count) / gamma(1)  ] (1/a = 1/1 = 1)
+		if (depth == d){
+				//sum = log_gammaratio(a, count);
+				sum = log_gammaratio(alpha_v_pa, count);				
+				//cerr << "count1: " << count << endl;
+
+		}
+		else if (depth <= d - 1){
+			for (int v = 0; v < arity; v ++){
+				if (children[v] != NULL){
+					double value = children[v]->evaluateHR(d, alpha_v_pa, r);
+					sum += value;
+				}
+			}
+			if (depth == d - 1){
+				int alpha_pa = alpha_v_pa / r;
+
+				sum -= log_gammaratio(alpha_pa, count);
+//				totalOf1 << "\nalpha_pa: " << alpha_pa << endl;
+//				totalOf1 << "count2: " << count << endl;
+//				totalOf1 << "-log_gammaratio(alpha_pa, count) = " << -log_gammaratio(alpha_pa, count) << endl;
+				//cerr << "count2: " << count << endl;
+				
+			}
+		}
+		return sum;
+	}
+	
+};
+
+
+class Model{
+public:
+	Model(){}
+	~Model(){}
+
+	void makeADTree(vector< vector<int> > & dm, vector<int> & arities){
+		//Create dmIndex from dm, dmIndex stroes {0, 1, R-1}
+		vector<int> dmIndex;
+//		cout << "dm.size() = " << dm.size() << endl;
+		for(int i = 0; i < (int) dm.size(); i++){
+			dmIndex.push_back(i);
+		}
+//		cout << "dmIndex:" << endl;
+//		for(int i = 0; i < (int) dm.size(); i++){
+//			cout << "dmIndex[" << i << "] = "<< dmIndex[i] << endl;
+//		}
+		
+		//  M: the num of attr
+//		int M = dm[0].size();
+//		cout << "M = " << M << endl;
+		
+		//cout << "Start MakeADTree()" << endl;
+		//ADTreeNode * ADTreeNodeP1 = MakeADTree(0, dmIndex, dm, arities);
+		ADTreeNodeP1 = MakeADTree(0, dmIndex, dm, arities);
+   		//cout << "End MakeADTree()" << endl;
+
+		//cout << "Start PrintADTreeNode()" << endl;
+		//PrintADTreeNode(ADTreeNodeP1, 0, M, arities);
+		//cout << "End PrintADTreeNode()" << endl;
+		
+	}
+
+	void makeADTree(){
+		makeADTree(data.dm, data.arities);
+	}
+
+	void freeADTree(vector< vector<int> > & dm, vector<int> & arities){
+		int M = dm[0].size();
+		//cout << "M = " << M << endl;
+		//cout << "Start FreeADTreeNode()" << endl;
+		FreeADTreeNode(ADTreeNodeP1, M, arities);
+		//cout << "End FreeADTreeNode()" << endl;
+		
+	}
+	
+	void freeADTree(){
+		freeADTree(data.dm, data.arities);
+	}
+	
+
+	void MakeADContiTable(vector< vector<int> > & dm, vector<int> & arities){
+		cout << "\nStart MakeADContiTable(vector< vector<int> > & dm, vector<int> & arities)" << endl;
+		//  M: the num of attr
+		int M = dm[0].size();
+		cout << "M = " << M << endl;
+		
+		deque<int> inquiryAttrs;
+	
+		for(int i = 0; i < (int) inquiryAttrs.size(); i++){
+			cout << "inquiryAttrs[" << i << "] = " << inquiryAttrs[i] << endl;
+		}
+		
+		vector<int> givenAttrs;
+		
+		vector<int> givenVals;
+			
+		cout << "Call MakeContab() " << endl;
+		CondiContiTable resultTable = MakeContab(inquiryAttrs, ADTreeNodeP1, givenAttrs, givenVals, arities);
+		cout << "Return MakeContab() " << endl;
+		for(unsigned int i = 0; i < resultTable.contiTableCounts.size(); i++){
+			cout << "resultTable.contiTableCounts[" << i << "] = " << resultTable.contiTableCounts[i] << endl; 
+		}
+		
+		cout << "End MakeADContiTable(vector< vector<int> > & dm, vector<int> & arities)" << endl;
+
+		
+	} //end void MakeADContiTable()
+	
+	void MakeADContiTable(){
+		MakeADContiTable(data.dm, data.arities);
+	}
+	
+	vector<int> MakeADContiTable(deque<int> inquiryAttrs, vector< vector<int> > & dm, vector<int> & arities){
+		//cout << "\nStart MakeADContiTable(vector< vector<int> > & dm, vector<int> & arities)" << endl;
+		//  M: the num of attr
+		//int M = dm[0].size();
+		//cout << "M = " << M << endl;
+		
+		//deque<int> inquiryAttrs;
+		
+		//Full inquiry
+//		for(int i = 0; i < M; i++){
+//			inquiryAttrs.push_back(i);
+//		}
+		
+		
+//		for(int i = 0; i < (int) inquiryAttrs.size(); i++){
+//			cout << "inquiryAttrs[" << i << "] = " << inquiryAttrs[i] << endl;
+//		}
+		
+		vector<int> givenAttrs;
+		
+		vector<int> givenVals;
+			
+		//cout << "Call MakeContab() " << endl;
+		CondiContiTable resultTable = MakeContab(inquiryAttrs, ADTreeNodeP1, givenAttrs, givenVals, arities);
+		//cout << "Return MakeContab() " << endl;
+//		
+//		for(unsigned int i = 0; i < resultTable.contiTableCounts.size(); i++){
+//			cout << "resultTable.contiTableCounts[" << i << "] = " << resultTable.contiTableCounts[i] << endl; 
+//		}
+//		
+//		cout << "End MakeADContiTable(vector< vector<int> > & dm, vector<int> & arities)" << endl;
+		return resultTable.contiTableCounts;
+		
+	} //end vector<int> MakeADContiTable()
+	
+
+	
+	void init(){
+		data.init();
+		m = data.numrecords;
+		n = data.numattributes;
+
+		layering.init(data);
+			
+		for (int a = 1; a < MAX_RQ; a ++){
+			double aa = 1.0 / a;
+			lr[a][0] = 0; lr2[a][0] = 0;
+			
+			for (int k = 1; k < MAX_COUNT; k ++){
+				lr[a][k] = lr[a][k-1] + log(aa + k - 1);
+				lr2[a][k] = lr2[a][k-1] + log(a + k - 1);
+			}
+		}
+	
+	}
+	
+	
+	
+	
+	
+	
+	
+	void test(){
+		for (int i = 0; i < 1; i ++){
+			vector<int> T;
+			T.clear();			
+			T.push_back(1);			
+			T.push_back(2);
+			
+		}
+		exit(1);
+	}
+	
+	
+	void testHR(){
+		for (int i = 0; i < 1; i ++){
+			vector<int> T;
+			T.clear();
+			T.push_back(1);
+			T.push_back(2);
+			T.push_back(3);
+			T.push_back(4);
+			
+		}
+	} //end void testHR()
+	
+	
+	void testHR_ADtree(){
+		for (int i = 0; i < 1; i ++){
+			deque<int> T;
+			T.clear();
+			T.push_back(1);
+			T.push_back(2);
+			T.push_back(3);
+			T.push_back(4);
+			
+		}
+	} //end void testHR_ADtree()
+	
+	
+
+	//--- Interface functions --- below ----------
+	//
+	// Log of local conditional probability, log p(xi | xS).
+	double log_lcp(int i, vector<int> & T){
+		
+		return log_lcp_mult(i, T);
+	}
+	double log_lcp(int i, int *T, int d){
+		vector<int> TT;
+		for (int j = 0; j < d; j ++) TT.push_back(T[j]);
+		return log_lcp_mult(i, TT);
+	}
+	
+
+	double log_lcpHR(int i, vector<int> & T){
+		
+		return log_lcp_multHR(i, T);
+	}
+	
+	double log_lcpHR(int i, int *T, int d){
+		vector<int> TT;
+		for (int j = 0; j < d; j ++) TT.push_back(T[j]);
+		//cerr << "start log_lcp_multHR(i, TT)" << endl;
+		double result = log_lcp_multHR(i, TT);
+		return result;
+	}
+	
+	
+	double log_lcpHR_ADtree(int i, int *T, int d){
+		deque<int> TT;
+		for (int j = 0; j < d; j ++) TT.push_back(T[j]);
+		double result = log_lcp_multHR_ADtree(i, TT);
+		return result;
+	}
+	
+	
+	double log_prior(int i, int *T, int d){
+		return lr[1][d] + lr[1][n-1-d] - lr[1][n-1];
+	}
+	int num_layers(){
+		return layering.numlayers;
+	}
+	void layer(int h, int** Vh, int* nh){
+		*Vh = &(layering.V[layering.cnh[h]]);
+		*nh = layering.nh[h];
+	}
+	void upper_layers(int h, int** Vu, int* nu){//similar to above but upper means parent
+		if (h < layering.numlayers){
+			*Vu = &(layering.V[layering.cnh[h+1]]);
+			*nu = n - layering.cnh[h+1];
+		}
+		else{ *Vu = NULL; *nu = 0;}	
+	}
+	void lower_layers(int h, int **Vl, int *nl){//similar to above but lower means children
+		if (h > 0){
+			*Vl = layering.V;
+			*nl = layering.cnh[h];
+		}
+		else{ *Vl = NULL; *nl = 0;}
+	}
+	int max_indegree(){
+		return atoi(Arguments::maxindegree);
+	}
+	bool edges_within(int h){
+		return layering.edges_within(h);
+	}
+	int num_nodes(){
+		return n;
+	}
+	
+	void print_edge_prob(ostream& f, int i, int j, double p){
+		//f<<data.heads[i]<<" -> "<<data.heads[j]<<" "<<p<<endl; 
+		f<<data.heads[i]<<" -> "<<data.heads[j]<<" \t" <<p <<endl; 
+	} 
+	
+	//
+	//--- Interface functions --- above ----------	
+	
+	
+private:	
+	// Dirichlet-multinomial model
+	double log_lcp_mult(int i, vector<int> & T){
+		vector<int> S;
+
+		S.push_back(i);
+				
+		for (int k = 0; k < (int)T.size(); k ++){
+			S.push_back(T[k]);
+		}
+		
+		//cerr << " Evaluation yields: " << evaluate(S) << " , log_gammaratio(1, 2) = " 
+		//	<< log_gammaratio(1, 2) << " , a = " <<  data.arities[S[0]] << endl;
+		
+		return evaluate(S);	
+	}
+	
+	double evaluate(vector<int> & S){
+		
+		Tnode root(data.maxarity+1, 0);
+		vector<int> u;
+		
+		//for each record t 
+		for (int t = 0; t < m; t ++){
+			u.clear();
+			for (int k = 0; k < (int)S.size(); k ++){
+				int j = S[k];
+				u.push_back(data.dm[t][j]);
+			}
+			//cerr << t+1 << "th update: ";
+			
+			root.update(S, u);
+			//cerr << endl;
+			
+		}
+		//root.print(cerr);
+		
+	
+		return root.evaluate(S.size(), 1, -data.arities[S[0]]);
+		
+	}
+	
+
+	double log_lcp_multHR(int i, vector<int> & T){
+		vector<int> S;
+				
+		S.push_back(i);
+				
+		for (int k = 0; k < (int)T.size(); k ++){
+			S.push_back(T[k]);
+		}
+	
+		//cerr << " Evaluation yields: " << evaluate(S) << " , log_gammaratio(1, 2) = " 
+		//	<< log_gammaratio(1, 2) << " , a = " <<  data.arities[S[0]] << endl;
+		double result = evaluateHR(S);
+		return result;
+		
+			
+	}//end log_lcp_multHR()
+	
+
+	double log_lcp_multHR_ADtree(int i, deque<int> & T){
+		//totalOf1 << "Inside log_lcp_multHR()" << endl;
+		
+		double result2 = evaluateHR_ADtree(T);
+
+		int insertIndex = -1;
+		for(int k = 0; k < (int)T.size(); k++){
+			if(i < T[k]){
+				insertIndex = k;
+				break;
+			}
+		}
+		if(insertIndex == -1){
+			insertIndex = (int)T.size();
+		}
+
+		deque<int> S;
+		
+		for(int k = 0; k <= insertIndex - 1; k++){
+			S.push_back(T[k]);
+		}
+		S.push_back(i);
+		for(int k = insertIndex; k < (int)T.size(); k++){
+			S.push_back(T[k]);
+		}
+		
+		double result1 = evaluateHR_ADtree(S);
+		
+		//cerr << " Evaluation yields: " << evaluate(S) << " , log_gammaratio(1, 2) = " 
+		//	<< log_gammaratio(1, 2) << " , a = " <<  data.arities[S[0]] << endl;
+		double result = result1 - result2;
+		return result;
+				
+	}
+	
+	
+	// Dirichlet-multinomial model
+	double evaluateHR(vector<int> & S){
+		Tnode * root = new Tnode (data.maxarity+1, 0);
+		vector<int> u;
+		
+		//for each record t 
+		for (int t = 0; t < m; t ++){
+			u.clear();
+			for (int k = 0; k < (int)S.size(); k ++){
+				int j = S[k];
+				u.push_back(data.dm[t][j]);
+			}
+			//cerr << t+1 << "th update: ";
+
+			(*root).update(S, u);
+			//
+			//cerr << endl;
+			
+		}
+		//root.print(cerr);
+		
+		int alpha_v_pa = 1;
+		for(int index = 0; index < (int)S.size(); index++){
+			alpha_v_pa *= data.arities[S[index]];
+		}
+		int r = data.arities[S[0]];
+		//return root.evaluateHR(S.size(), alpha_v_pa, r);
+
+		double result = (*root).evaluateHR(S.size(), alpha_v_pa, r);
+		//totalOf1 << "    result of root.evaluateHR(" << S.size() <<", " 
+		//	<< alpha_v_pa << ", " << r << ") = " << result << endl;
+		delete root;
+		return result;
+		
+	} //end double evaluateHR(vector<int> & S)
+	
+
+	double evaluateHR_ADtree(deque<int> & S){
+
+		int alpha = 1;
+		for(int index = 0; index < (int)S.size(); index++){
+			alpha *= data.arities[S[index]];
+		}
+		
+		vector<int> countVec = MakeADContiTable(S, data.dm, data.arities);
+		
+		double sum = 0;
+		for(int i = 0; i < (int) countVec.size(); i++){
+			//totalOf1 << "countVec[" << i <<" ]" << countVec[i] << endl;
+			sum += log_gammaratio(alpha, countVec[i]);
+		}
+		
+//		totalOf1 << "evaluateHR_ADtree(deque<int> & S) ends" << endl;
+
+		return sum;
+		
+	} //end double evaluateHR_ADtree()
+	
+		
+	Data data;
+	int n;
+	int m;
+	Layering layering;	
+	
+	ADTreeNode * ADTreeNodeP1;
+	
+}; //end of class Model
+
+
+
+
+
+#endif