//---------------------------------------------------------------------------

#include "stdafx.h"

#ifdef _MSC_VER
#pragma hdrstop
#endif

#include "kernelincludes.h"

#include "orthogonalize.h"
//---------------------------------------------------------------------------

Ortho::Ortho( PLA *_pla )
{
	pla = _pla;
    ok_terms = new PLA(pla);
	if ( ok_terms == NULL) throw RuntimeException( "Ortho: Not enough memory" );
	ok_terms->SetUnique();
}

void Ortho::Orthogonalize(void)
/* Randomly select one of two intersecting terms, randomly splits one
*/
{
    string s, s1;

    Term newterm;
    list<Term>::iterator ti1, ti2, ti3, ti4;
    vector<int> candidates;
	vector<int>::iterator ci;
    unsigned int pos;
    unsigned int i, j;
    list<Term> *data;

	data = pla->GetData();
			// first preprocessing
	RemoveAbsorbed();
	pla->MergeTerms();
	MoveOKTerms2();

    ti1 = data->begin();
    while ( ti1 != data->end()) {
       ti2 = ti1;
	   ti2++;
	   do {
            if ( ti2 == data->end() ) break;
			j = (*ti1).Relation2(*ti2);
			if ( j == 0 || (j & 8) == 0 ) ti2++;				// no IM or OM intersection
			else if ( j == 63 || j == 27 )
				ti2 = pla->DeleteTerm(ti2);		// equal or contained
			else if ( (j & 7) == 7 ) {				// equal inputs, different outputs. Merge them
                (*ti1).Merge(*ti2);
                ti2 = pla->DeleteTerm(ti2);
			} else if ( (j & 7) != 0 && (j & 8) != 0 ) break;   // non-empty IM and OM intersection
			else throw RuntimeException("Orthogonalize: Something's wrong");
       } while ( ti2 != data->end() );
       if ( ti2 != data->end()) {       // they do intersect
            for ( i = 0; i < pla->GetInputs(); i++ )
                if ( (*ti1).GetIM(i) == '-' && (*ti2).GetIM(i) != '-' ) candidates.push_back(i + 1);
                else if ( (*ti2).GetIM(i) == '-' && (*ti1).GetIM(i) != '-' ) candidates.push_back(-(int)(i + 1));
            if ( candidates.size() == 0 ) throw RuntimeException("Orthogonalize: Something's wrong");
            else {
				i = rand() % candidates.size();
				if ( candidates[i] > 0 ) ti3 = ti1; else ti3 = ti2;
                pos = abs(candidates[i]) - 1;
                candidates.clear();
                newterm = *ti3;         // make a copy
                (*ti3).SetIM(pos, '0');   // set 0 in the original term
                ti4 = data->begin();
				do {
					if ( ti4 != ti3 ) { 
                        if ( (*ti4).Contains(*ti3) ) {
                            pla->DeleteTerm(ti3);
							ti1 = data->begin();
                            break;
                        } else if ( (*ti3).Contains(*ti4) ) {
                            ti4 = pla->DeleteTerm(ti4);
							ti1 = data->begin();
                        } else ti4++;                    
					} else ti4++;
				} while ( ti4 != data->end() );
                newterm.SetIM(pos, '1');    // set 1 in the new term
                for ( ti4 = data->begin(); ti4 != data->end(); ti4++ ) {
                        if ( (*ti4).Contains(newterm) ) break;	// already included somewhere
                }
                if ( ti4 == data->end() ) pla->AppendTerm(newterm);
            }
         //   pla->SavePLA("tmp.pla");
       } else ti1++;
    }
	MoveOKTerms2();
	if ( pla->GetTerms() != 0 ) throw RuntimeException( "Ortho::Orthogonalize:: Something's wrong - resulting PLA is not orthogonal" );
}

void Ortho::MoveOKTerms(void)
// Moves terms intersecting with no other terms to ok_terms. Matrix must be computed
{
	unsigned int i, j;
	list<Term>::iterator ti;

	for ( i = 0; i < pla->GetTerms(); i++ ) {
		for ( j = 0; j < pla->GetTerms(); j++ )
			if ( matrix[i][j] != -1 ) break;
		if ( j == pla->GetTerms() ) {
			ti = pla->GetTermByIndex(i);
			ok_terms->AppendTerm(*ti);
			(*ti).SetUnwanted();
		}
	}
	pla->Purge();
}

void Ortho::MoveOKTerms2(void)
// Moves terms intersecting with no other terms to ok_terms. Does not rely on matrix.
{
	list<Term>::iterator ti1, ti2;
	bool flag;

	for ( ti1 = pla->GetData()->begin(); ti1 != pla->GetData()->end(); ti1++ ) {
		flag = true;
		for ( ti2 = pla->GetData()->begin(); ti2 != pla->GetData()->end(); ti2++ ) 
			if ( ti1 != ti2 )
				if ( (*ti1).Intersects(*ti2)) {
					flag = false;
					break;
				}
		if ( flag ) {
			ok_terms->AppendTerm(*ti1);
			(*ti1).SetUnwanted();
		}
	}
	pla->Purge();
}

void Ortho::GenerateMatrix(void)
{
	vector<int> matrixrow;
    unsigned int i, j, k;
    list<Term>::iterator ti1, ti2;

	matrix.clear();
	matrixrow.clear();
	matrixrow.resize(pla->GetTerms(), 0);	
	for ( i = 0; i < pla->GetTerms(); i++ ) 
		matrix.push_back(matrixrow);
	for ( i = 0; i < pla->GetTerms(); i++ ) {
		matrix[i][i] = -1;
		ti1 = pla->GetTermByIndex(i);
		for ( j = i+1; j < pla->GetTerms(); j++ ) {
			ti2 = pla->GetTermByIndex(j);
			k = (*ti1).IntersectSize(*ti2);
			matrix[i][j] = k;
			matrix[j][i] = k;
		}
	}
}

void Ortho::RecomputeMatrix(void)
{
    unsigned int i, j, k;
    list<Term>::iterator ti1, ti2;

	for ( i = 0; i < pla->GetTerms(); i++ ) {
		matrix[i][i] = -1;
		ti1 = pla->GetTermByIndex(i);
		for ( j = i+1; j < pla->GetTerms(); j++ ) {
			ti2 = pla->GetTermByIndex(j);
			k = (*ti1).IntersectSize(*ti2);
			matrix[i][j] = k;
			matrix[j][i] = k;
		}
	}
}

void Ortho::AppendMatrixRow( const TermList::iterator ti )
/*
	Increases matrix by 1 and computes values for added term ti
*/
{
    vector<int> matrixrow;
	unsigned int i, j, k, l;
    list<Term>::iterator ti2;

	i = matrix.size();
	l = pla->GetIndex(ti);
	if ( i + 1 != pla->GetTerms() ) throw RuntimeException( "Ortho::AppendMatrixRow: Wrong matrix size" );
	matrixrow.resize(pla->GetTerms());	
	matrix.push_back(matrixrow);
	for ( j = 0; j < i; j++ ) {
		matrix[j].push_back(-1);
		ti2 = pla->GetTermByIndex(j);
		k = (*ti).IntersectSize(*ti2);
		matrix[i][j] = k;
		matrix[j][i] = k;
	}
}

void Ortho::RecomputeMatrixRow( const TermList::iterator ti )
/*
	Recomputes values for term ti
*/
{
    unsigned int i, j, k;
    list<Term>::iterator ti2;

//	pla->GetTerm(0);
	i = pla->GetIndex(ti);
	for ( j = 0; j < pla->GetTerms(); j++ ) 
		if ( j != i && matrix[i][j] != -1 ) {
			ti2 = pla->GetTermByIndex(j);
			k = (*ti).IntersectSize(*ti2);
			matrix[i][j] = k;
			matrix[j][i] = k;
		}
}

void Ortho::RemoveAbsorbed(void)
//	Removed absorbed terms from PLA
//  If IM is equal to some existing term, merges the term and deletes ti.
{
    TermList::iterator ti, ti1;

	if ( pla->GetData()->size() < 2 ) return;
	ti = pla->GetData()->begin();
	do {
        ti1 = ti;
        ti1++;
        while ( ti1 != pla->GetData()->end() ) {
			if ( (*ti).Contains( *ti1 ) ) ti1 = pla->DeleteTerm( ti1 );
			else if ( (*ti1).Contains( *ti ) ) {
				ti = pla->DeleteTerm( ti );
				break;
			}
			else if ( (*ti).GetIM() == (*ti1).GetIM() ) {
				(*ti).Merge(*ti1);
				ti1 = pla->DeleteTerm( ti1 );
			} else ti1++;
        }
		ti++;
    } while ( ti != pla->GetData()->end() );
}

bool Ortho::RemoveAbsorbed( const TermList::iterator ti )
//	Remove ti if absorbed
//  If IM is equal to some existing term, merges the term and deletes ti.
//  Returns true if it was removed
{
    TermList::iterator ti1;

    for ( ti1 = pla->GetData()->begin(); ti1 != pla->GetData()->end(); ti1++ ) 
		if ( ti1 != ti ) {
			if ( (*ti1).Contains( *ti ) ) {
				pla->DeleteTerm( ti );
				return true;
			}
			if ( (*ti).GetIM() == (*ti1).GetIM() ) {
				(*ti1).Merge(*ti);
				pla->DeleteTerm( ti );
				return true;
			}
		}
	return false;
}

bool Ortho::TryToMerge( const TermList::iterator ti )
//	Tries to merge the term
//  Returns true if it was merged
{
    TermList::iterator ti1;

    for ( ti1 = pla->GetData()->begin(); ti1 != pla->GetData()->end(); ti1++ ) 
		if ( ti1 != ti ) {
			if ( (*ti).GetIM() == (*ti1).GetIM() ) {
				(*ti1).Merge(*ti);
				return true;
			}
		}
	return false;
}

bool Ortho::IsAbsorbed( const Term &t ) const
//	Returns true if t is absorbed
//  If IM is equal to some existing term, merges the term and returns true.
{
    TermList::iterator ti1;

	for ( ti1 = pla->GetData()->begin(); ti1 != pla->GetData()->end(); ti1++ ) {
		if ( t.GetIM() == (*ti1).GetIM() ) {
			(*ti1).Merge(t);
			return true;
		}
		else if ( (*ti1).Contains(t) ) {
			return true;
		}
	}
	return false;
}

bool Ortho::WillBeAbsorbed( const TermList::iterator ti ) const
//	Returns true if t can be absorbed or merged
{
    TermList::iterator ti1;

	for ( ti1 = pla->GetData()->begin(); ti1 != pla->GetData()->end(); ti1++ ) 
		if ( ti != ti1 ) {
			if ( (*ti).GetIM() == (*ti1).GetIM() ) return true;
			else if ( (*ti1).Contains((*ti)) ) 
				return true;
			
		}
	return false;
}

void Ortho::Orthogonalize2(void)
/*  Select term intersecting with most terms.
	Then select splitting variables, so that one term (both) will be absorbed.
	From it, splits by variable present in most intersected terms
*/
{
	Term newterm;
    list<Term>::iterator ti1, ti2, ti3, ti4;
    vector<int> candidates;
	vector<splitvar> candidates2;
	vector<int> matrixrow;
    unsigned int prom;
	short pol;
    unsigned int i, j;
    list<Term> *data;
	unsigned int cnt;
	unsigned int mx, cnt1, cnt0;
	unsigned int term;
	unsigned int val;
	bool absorbed1, absorbed2;		// are old/new terms absorbed?
	bool abs2;					// both terms are absorbed

	unsigned int tr1, tr2, tr3;
	tr1 = 0;
	tr2 = 0;
	tr3 = 0;

//	ofstream debugfile;
//	debugfile.open("ortho_debug.txt", std::ios::app);

	data = pla->GetData();
				// first preprocessing
	RemoveAbsorbed();
	pla->MergeTerms();
	MoveOKTerms2();
//	pla->SavePLA("bad.pla");
//	ok_terms->SavePLA("ok.pla");
	GenerateMatrix();

	while ( pla->GetTerms() != 0 ) {
		
			// find term conflicting with most terms
		mx = 0;
		candidates.clear();
		for ( i = 0; i < pla->GetTerms(); i++ ) 
			if ( (*pla->GetTermByIndex(i)).GetDim() != 0 ) {
				cnt = 0;
				for ( j = 0; j < pla->GetTerms(); j++ ) if ( matrix[i][j] != -1 ) cnt++;
				if ( cnt > mx ) {
					candidates.clear();
					candidates.push_back(i);
					mx = cnt;
				} else if ( cnt == mx && cnt != 0 )
					candidates.push_back(i);
			}
				// select one randomly
		if ( candidates.size() == 0 ) throw RuntimeException( "Ortho::Something's wrong - no intersection terms found" );
		term = candidates[rand() % candidates.size()];
		ti1 = pla->GetTermByIndex(term);

				// select splitting variable, so that one term (both) will be absorbed
		abs2 = false;
		candidates2.clear();
		for ( i = 0; i < pla->GetInputs(); i++ ) 
			if ( (*ti1).GetIM(i) == '-' ) {
				(*ti1).SetIM( i, '1' );
				if ( abs2 ) {	// now both have to be absorbed
					if ( WillBeAbsorbed(ti1) ) {
						(*ti1).SetIM( i, '0' );
						if ( WillBeAbsorbed(ti1) ) {
							s_var.var = i;
							s_var.pol = 0;
							candidates2.push_back(s_var);
						}
					}
				} else {
					absorbed1 = WillBeAbsorbed(ti1);
					(*ti1).SetIM( i, '0' );
					absorbed2 = WillBeAbsorbed(ti1);
					if ( absorbed1 && absorbed2 ) {
						abs2 = 1;
						candidates.clear();
						s_var.var = i;
						s_var.pol = 0;
						candidates2.push_back(s_var);
					} else if ( absorbed1 ) {
						s_var.var = i;
						s_var.pol = 1;
						candidates2.push_back(s_var);
					} else if ( absorbed2 ) {
						s_var.var = i;
						s_var.pol = -1;
						candidates2.push_back(s_var);
					}
				}
				(*ti1).SetIM( i, '-' );
			}

			// if there are no candidates, select the variable present in most intersected terms
//		debugfile << candidates2.size() << endl;
		if ( candidates2.empty() ) {
			mx = 0;
			for ( i = 0; i < pla->GetInputs(); i++ ) {
				if ( (*ti1).GetIM(i) == '-' ) {
					cnt0 = 0;
					cnt1 = 0;
					for ( ti2 = data->begin(); ti2 != data->end(); ti2++ ) {
						if ( matrix[term][pla->GetIndex(ti2)] != -1 )
							if ( (*ti2).GetIM(i) == '0' ) cnt1++;
							else if ( (*ti2).GetIM(i) == '1' ) cnt0++;
					}
					val = max( cnt0, cnt1 );
					if ( mx < val ) {
						candidates2.clear();
						s_var.var = i;
						s_var.pol = 0;
						candidates2.push_back(s_var);
						mx = val;
					} else if ( mx == val ) {
						s_var.var = i;
						s_var.pol = 0;
						candidates2.push_back(s_var);
					}
				}
			}
		}
				// select one randomly
		if ( candidates2.empty() ) throw RuntimeException( "Ortho::Something's wrong - no intersection terms found" );
		i = rand() % candidates2.size();
		prom = candidates2[i].var;
		pol = candidates2[i].pol;
		
				// split the term
		newterm = *ti1;
		(*ti1).SetIM( prom, '1' );
		newterm.SetIM( prom, '0' );
		absorbed1 = RemoveAbsorbed( ti1 );
		absorbed2 = IsAbsorbed( newterm );
		if ( !absorbed2 ) ti2 = pla->AppendTerm( newterm );

			// recompute everything
		tr1++;
		if ( !absorbed1 && absorbed2 ) {
			RecomputeMatrix();
			tr2++;
		}
		// else if ( !absorbed1 && !absorbed2 ) AppendMatrixRow( ti2 );   // this only slows it down
		else {
			GenerateMatrix();
			tr3++;
		}
		i = pla->GetTerms();
		MoveOKTerms();
		if ( pla->GetTerms() != i ) GenerateMatrix();
//		pla->SavePLA("bad.pla");
//		ok_terms->SavePLA("ok.pla");
	}
	if ( !IsOrthogonal(ok_terms) ) throw RuntimeException( "Ortho::Something's wrong - resulting PLA is not orthogonal" );
//	debugfile << endl;
//	debugfile.close();
//	cout << "   All/Recomputed/Generated " << tr1 << " " << tr2 << " " << tr3 << endl;
}

void Ortho::Orthogonalize2a(void)
/*  Select term intersecting with most terms.
	Then select splitting variables, so that one term (both) will be absorbed.
	From it, splits by variable present in most intersected terms
	No intersection matrix is used.

	It is slower than Orthogonalize2!
*/
{
	Term newterm;
    list<Term>::iterator ti1, ti2, ti3, ti4;
    vector< list<Term>::iterator > candidates;	// candidate terms for the selection
	vector<splitvar> candidates2;	// candidate variables for the selection
    unsigned int prom;
	short pol;
    unsigned int i;
    list<Term> *data;
	unsigned int cnt;
	unsigned int mx, cnt1, cnt0;
	unsigned int val;
	bool absorbed1, absorbed2;		// are old/new terms absorbed?
	bool abs2;					// both terms are absorbed

//	ofstream debugfile;
//	debugfile.open("ortho_debug.txt", std::ios::app);

	data = pla->GetData();
				// first preprocessing
	RemoveAbsorbed();
	pla->MergeTerms();
	MoveOKTerms2();
//	pla->SavePLA("bad.pla");
//	ok_terms->SavePLA("ok.pla");

	 while ( pla->GetTerms() != 0 ) {
			// find term conflicting with most terms
		mx = 0;
		candidates.clear();
		for ( ti1 = pla->GetData()->begin(); ti1 != pla->GetData()->end(); ti1++ ) 
			if ( (*ti1).GetDim() != 0 ) {
				cnt = 0;
				ti2 = ti1;
				ti2++;
				for ( ti2 = pla->GetData()->begin(); ti2 != pla->GetData()->end(); ti2++ )
					if ( ti1 != ti2 ) 
						if ( (*ti1).Intersects(*ti2) ) cnt++;
				if ( cnt > mx ) {
					candidates.clear();
					candidates.push_back(ti1);
					mx = cnt;
				} else if ( cnt == mx && cnt != 0 )
					candidates.push_back(ti1);
			}
				// select one randomly
		if ( candidates.size() == 0 ) throw RuntimeException( "Ortho::Something's wrong - no intersection terms found" );
		ti1 = candidates[rand() % candidates.size()];

				// select splitting variable, so that one term (both) will be absorbed
		abs2 = false;
		candidates2.clear();
		for ( i = 0; i < pla->GetInputs(); i++ ) 
			if ( (*ti1).GetIM(i) == '-' ) {
				(*ti1).SetIM( i, '1' );
				if ( abs2 ) {	// now both have to be absorbed
					if ( WillBeAbsorbed(ti1) ) {
						(*ti1).SetIM( i, '0' );
						if ( WillBeAbsorbed(ti1) ) {
							s_var.var = i;
							s_var.pol = 0;
							candidates2.push_back(s_var);
						}
					}
				} else {
					absorbed1 = WillBeAbsorbed(ti1);
					(*ti1).SetIM( i, '0' );
					absorbed2 = WillBeAbsorbed(ti1);
					if ( absorbed1 && absorbed2 ) {
						abs2 = 1;
						candidates.clear();
						s_var.var = i;
						s_var.pol = 0;
						candidates2.push_back(s_var);
					} else if ( absorbed1 ) {
						s_var.var = i;
						s_var.pol = 1;
						candidates2.push_back(s_var);
					} else if ( absorbed2 ) {
						s_var.var = i;
						s_var.pol = -1;
						candidates2.push_back(s_var);
					}
				}
				(*ti1).SetIM( i, '-' );
			}

			// if there are no candidates, select the variable present in most intersected terms
//		debugfile << candidates2.size() << endl;
		if ( candidates2.empty() ) {
			mx = 0;
			for ( i = 0; i < pla->GetInputs(); i++ ) {
				if ( (*ti1).GetIM(i) == '-' ) {
					cnt0 = 0;
					cnt1 = 0;
					for ( ti2 = data->begin(); ti2 != data->end(); ti2++ ) {
						if ( (*ti1).Intersects(*ti2) )
							if ( (*ti2).GetIM(i) == '0' ) cnt1++;
							else if ( (*ti2).GetIM(i) == '1' ) cnt0++;
					}
					val = max( cnt0, cnt1 );
					if ( mx < val ) {
						candidates2.clear();
						s_var.var = i;
						s_var.pol = 0;
						candidates2.push_back(s_var);
						mx = val;
					} else if ( mx == val ) {
						s_var.var = i;
						s_var.pol = 0;
						candidates2.push_back(s_var);
					}
				}
			}
		}
				// select one randomly
		if ( candidates2.empty() ) throw RuntimeException( "Ortho::Something's wrong - no intersection terms found" );
		i = rand() % candidates2.size();
		prom = candidates2[i].var;
		pol = candidates2[i].pol;
		
				// split the term
		newterm = *ti1;
		(*ti1).SetIM( prom, '1' );
		newterm.SetIM( prom, '0' );
		absorbed1 = RemoveAbsorbed( ti1 );
		absorbed2 = IsAbsorbed( newterm );
		if ( !absorbed2 ) ti2 = pla->AppendTerm( newterm );

			// recompute everything
		MoveOKTerms2();
//		pla->SavePLA("bad.pla");
//		ok_terms->SavePLA("ok.pla");
	}
//	if ( !IsOrthogonal(ok_terms) ) throw RuntimeException( "Ortho::Something's wrong - resulting PLA is not orthogonal" );
//	debugfile << endl;
//	debugfile.close();
}

void Ortho::Orthogonalize3(void)
/*  1. find a term & variable to split, so that one (or both) split terms will be absorbed
	2. from these candidates (if any) find term conflicting with most terms
	3. from these candidates select the variable present in most intersected terms
*/
{
	Term newterm;
    list<Term>::iterator ti1, ti2, ti3, ti4;
    vector<int> candidates;
	vector<splitvar> candidates2;
	vector<splitterm> candidates3;
	vector<splitterm>::iterator cand3_it;
	vector<int> matrixrow;
    unsigned int prom;
	short pol;
    unsigned int i, j;
    list<Term> *data;
	unsigned int cnt;
	unsigned int mx, cnt1, cnt0;
	unsigned int term;
	unsigned int val;
	bool absorbed1, absorbed2;		// are old/new terms absorbed?
	bool abs2;					// both terms are absorbed
	set<unsigned int> termset;
	set<unsigned int> varset;
	set<unsigned int>::iterator varset_it;

	unsigned int tr1, tr2, tr3;
	tr1 = 0;
	tr2 = 0;
	tr3 = 0;

//	ofstream debugfile;
//	debugfile.open("ortho_debug.txt", std::ios::app);

	data = pla->GetData();
				// first preprocessing
	RemoveAbsorbed();
	pla->MergeTerms();
	MoveOKTerms2();
//	pla->SavePLA("bad.pla");
//	ok_terms->SavePLA("ok.pla");
	GenerateMatrix();

	while ( pla->GetTerms() != 0 ) {

			// find a term & variable to split, so that one (or both) split terms will be absorbed
		abs2 = false;
		candidates3.clear();
		termset.clear();
		for ( ti1 = pla->GetData()->begin(); ti1 != pla->GetData()->end(); ti1++ ) 
			if ( (*ti1).GetDim() != 0 )
				for ( i = 0; i < pla->GetInputs(); i++ ) 
					if ( (*ti1).GetIM(i) == '-' ) {
						(*ti1).SetIM( i, '1' );
						if ( abs2 ) {	// now both have to be absorbed to be accepted
							if ( WillBeAbsorbed(ti1) ) {
								(*ti1).SetIM( i, '0' );
								if ( WillBeAbsorbed(ti1) ) {
									s_term.ti = ti1;
									s_term.var = i;
									s_term.pol = 0;
									candidates3.push_back(s_term);
									termset.insert(pla->GetIndex(ti1));
								}
							}
						} else {
							absorbed1 = WillBeAbsorbed(ti1);
							(*ti1).SetIM( i, '0' );
							absorbed2 = WillBeAbsorbed(ti1);
							if ( absorbed1 && absorbed2 ) {
								abs2 = true;
								candidates3.clear();
								s_term.ti = ti1;
								s_term.var = i;
								s_term.pol = 0;
								termset.insert(pla->GetIndex(ti1));
								candidates3.push_back(s_term);
							} else if ( absorbed1 ) {
								s_term.ti = ti1;
								s_term.var = i;
								s_term.pol = 1;
								termset.insert(pla->GetIndex(ti1));
								candidates3.push_back(s_term);
							} else if ( absorbed2 ) {
								s_term.ti = ti1;
								s_term.var = i;
								s_term.pol = -1;
								termset.insert(pla->GetIndex(ti1));
								candidates3.push_back(s_term);
							}
						}
						(*ti1).SetIM( i, '-' );
					}

					// from these candidates (if any) find term conflicting with most terms
		mx = 0;
		candidates.clear();
		for ( i = 0; i < pla->GetTerms(); i++ ) 
			if ( (*pla->GetTermByIndex(i)).GetDim() != 0 && (termset.empty() || termset.find(i) != termset.end()) ) {
				cnt = 0;
				for ( j = 0; j < pla->GetTerms(); j++ ) if ( matrix[i][j] != -1 ) cnt++;
				if ( cnt > mx ) {
					candidates.clear();
					candidates.push_back(i);
					mx = cnt;
				} else if ( cnt == mx && cnt != 0 )
					candidates.push_back(i);
			}
				// select one randomly
		if ( candidates.size() == 0 ) throw RuntimeException( "Orthogonalize3::Something's wrong - no intersection terms found" );
		term = candidates[rand() % candidates.size()];
		ti1 = pla->GetTermByIndex(term);

			// create a candidate set for splitting variables
		varset.clear();
		if ( candidates3.empty() ) {
			for ( i = 0; i < pla->GetInputs(); i++ ) 
				if ( (*ti1).GetIM(i) == '-' ) varset.insert(i);
		} else {
			for ( cand3_it = candidates3.begin(); cand3_it != candidates3.end(); cand3_it++ )
				if ( (*cand3_it).ti == ti1 ) varset.insert( (*cand3_it).var );
		}

			// select the variable present in most intersected terms
		if ( varset.size() > 1 ) {
			mx = 0;
			candidates2.clear();
			for ( varset_it = varset.begin();  varset_it != varset.end(); varset_it++ ) {
				i = *varset_it;
				cnt0 = 0;
				cnt1 = 0;
				for ( ti2 = data->begin(); ti2 != data->end(); ti2++ ) {
					if ( matrix[term][pla->GetIndex(ti2)] != -1 )
						if ( (*ti2).GetIM(i) == '0' ) cnt1++;
						else if ( (*ti2).GetIM(i) == '1' ) cnt0++;
				}
				val = max( cnt0, cnt1 );
				if ( mx < val ) {
					candidates2.clear();
					s_var.var = i;
					s_var.pol = 0;
					candidates2.push_back(s_var);
					mx = val;
				} else if ( mx == val ) {
					s_var.var = i;
					s_var.pol = 0;
					candidates2.push_back(s_var);
				}
			}
					// select one randomly
//			pla->SavePLA("bad.pla");
			if ( candidates2.empty() ) throw RuntimeException( "Orthogonalize3::Something's wrong - no intersection terms found" );
			i = rand() % candidates2.size();
			prom = candidates2[i].var;
			pol = candidates2[i].pol;
		} else if (varset.size() == 1 ) {
			prom = *(varset.begin());
			pol = 0;
		} else {						// varset empty. No splitting variable helps. Select one randomly
			i = rand() % pla->GetInputs() + 1;
			prom = 0;
			do {
				if ( (*ti1).GetIM(prom) == '-' ) i--;
				else {
					prom++;
					prom %= pla->GetInputs();
				}
			} while ( i > 0 );
		}


//		debugfile << candidates3.size() << endl;
		
				// split the term
		if ( (*ti1).GetIM(prom) != '-' ) throw RuntimeException( "Orthogonalize3::Something's wrong - bad splitting variable" );
		newterm = *ti1;
		(*ti1).SetIM( prom, '1' );
		newterm.SetIM( prom, '0' );
		absorbed1 = RemoveAbsorbed( ti1 );
		absorbed2 = IsAbsorbed( newterm );
		if ( !absorbed2 ) ti2 = pla->AppendTerm( newterm );

			// recompute matrix
		tr1++;

		if ( !absorbed1 && absorbed2 ) {
			RecomputeMatrix();
			tr2++;
		}
		// else if ( !absorbed1 && !absorbed2 ) AppendMatrixRow( ti2 );   // this only slows it down
		else {
			GenerateMatrix();
			tr3++;
		}
		i = pla->GetTerms();
		MoveOKTerms();
		if ( pla->GetTerms() != i && pla->GetTerms() != 0 ) GenerateMatrix();
		// pla->SavePLA("bad.pla");
		// ok_terms->SavePLA("ok.pla");
	}
	if ( !IsOrthogonal(ok_terms) ) throw RuntimeException( "Orthogonalize3::Something's wrong - resulting PLA is not orthogonal" );
//	debugfile << endl;
//	debugfile.close();
//	cout << "   All/Recomputed/Generated: " << tr1 << " " << tr2 << " " << tr3 << endl;
}

void Ortho::Orthogonalize4(void)
/*  1. find a term & variable to split, so that one (or both) split terms will be absorbed.
	2. from these, choose terms&variables, where splitting variables are present in most intersected terms
*/
{
	Term newterm;
    list<Term>::iterator ti1, ti2, ti3, ti4;
    vector<int> candidates;
	vector<splitterm> candidates3;
	vector<splitterm> candidates2;
	vector<splitterm>::iterator cand3_it;
	vector<int> matrixrow;
    unsigned int prom;
    unsigned int i;
    list<Term> *data;
	unsigned int mx, cnt1, cnt0;
	unsigned int val;
	bool absorbed1, absorbed2;		// are old/new terms absorbed?
	bool abs2;					// both terms are absorbed

	unsigned int tr1, tr2, tr3;
	tr1 = 0;
	tr2 = 0;
	tr3 = 0;

//	ofstream debugfile;
//	debugfile.open("ortho_debug.txt", std::ios::app);

	data = pla->GetData();
				// first preprocessing
	RemoveAbsorbed();
	pla->MergeTerms();
	MoveOKTerms2();
//	pla->SavePLA("bad.pla");
//	ok_terms->SavePLA("ok.pla");
	GenerateMatrix();

	while ( pla->GetTerms() != 0 ) {
			// find a term & variable to split, so that one (or both) split terms will be absorbed
		abs2 = false;
		candidates3.clear();
		for ( ti1 = pla->GetData()->begin(); ti1 != pla->GetData()->end(); ti1++ ) 
			if ( (*ti1).GetDim() != 0 )
				for ( i = 0; i < pla->GetInputs(); i++ ) 
					if ( (*ti1).GetIM(i) == '-' ) {
						(*ti1).SetIM( i, '1' );
						if ( abs2 ) {	// now both have to be absorbed to be accepted
							if ( WillBeAbsorbed(ti1) ) {
								(*ti1).SetIM( i, '0' );
								if ( WillBeAbsorbed(ti1) ) {
									s_term.ti = ti1;
									s_term.var = i;
									s_term.pol = 0;
									candidates3.push_back(s_term);
								}
							}
						} else {
							absorbed1 = WillBeAbsorbed(ti1);
							(*ti1).SetIM( i, '0' );
							absorbed2 = WillBeAbsorbed(ti1);
							if ( absorbed1 && absorbed2 ) {
								abs2 = true;
								candidates3.clear();
								s_term.ti = ti1;
								s_term.var = i;
								s_term.pol = 0;
								candidates3.push_back(s_term);
							} else if ( absorbed1 ) {
								s_term.ti = ti1;
								s_term.var = i;
								s_term.pol = 1;
								candidates3.push_back(s_term);
							} else if ( absorbed2 ) {
								s_term.ti = ti1;
								s_term.var = i;
								s_term.pol = -1;
								candidates3.push_back(s_term);
							}
						}
						(*ti1).SetIM( i, '-' );
					}
		if ( candidates3.empty() ) 		// no candidates, generate all (terms x DC-vars)	
			for ( ti1 = pla->GetData()->begin(); ti1 != pla->GetData()->end(); ti1++ ) 
				if ( (*ti1).GetDim() != 0 )
					for ( i = 0; i < pla->GetInputs(); i++ )
						if ( (*ti1).GetIM(i) == '-' ) {
							s_term.ti = ti1;
							s_term.var = i;
							s_term.pol = 0;
							candidates3.push_back(s_term);
						}
		if ( candidates3.empty() ) throw RuntimeException( "Orthogonalize::Ortho4: Something's wrong. No splitting variables" );

		// select the variable (+term) present in most intersected terms
			
		mx = 0;
		candidates2.clear();
		for ( cand3_it = candidates3.begin();  cand3_it != candidates3.end(); cand3_it++ ) {
			cnt0 = 0;
			cnt1 = 0;
			for ( ti2 = data->begin(); ti2 != data->end(); ti2++ ) {
				if ( matrix[pla->GetIndex((*cand3_it).ti)][pla->GetIndex(ti2)] != -1 )
					if ( (*ti2).GetIM((*cand3_it).var) == '0' ) cnt1++;
					else if ( (*ti2).GetIM((*cand3_it).var) == '1' ) cnt0++;
			}
			val = max( cnt0, cnt1 );
			if ( mx < val ) {
				candidates2.clear();
				candidates2.push_back(*cand3_it);
				mx = val;
			} else if ( mx == val && mx != 0 )
				candidates2.push_back(*cand3_it);
		}
		if ( candidates2.empty() ) {		// splitting these terms by its variables will not help. Do it once more with all terms.
			for ( ti1 = pla->GetData()->begin(); ti1 != pla->GetData()->end(); ti1++ ) 
				for ( i = 0; i < pla->GetInputs(); i++ )
					if ( (*ti1).GetIM(i) == '-' ) {
						s_term.ti = ti1;
						s_term.var = i;
						s_term.pol = 0;
						candidates3.push_back(s_term);
					}
			for ( cand3_it = candidates3.begin();  cand3_it != candidates3.end(); cand3_it++ ) {
				cnt0 = 0;
				cnt1 = 0;
				for ( ti2 = data->begin(); ti2 != data->end(); ti2++ ) {
					if ( matrix[pla->GetIndex((*cand3_it).ti)][pla->GetIndex(ti2)] != -1 )
						if ( (*ti2).GetIM((*cand3_it).var) == '0' ) cnt1++;
						else if ( (*ti2).GetIM((*cand3_it).var) == '1' ) cnt0++;
				}
				val = max( cnt0, cnt1 );
				if ( mx < val ) {
					candidates2.clear();
					candidates2.push_back(*cand3_it);
					mx = val;
				} else if ( mx == val && mx != 0 )
					candidates2.push_back(*cand3_it);
			}
		}
		if ( candidates2.empty() ) throw RuntimeException( "Orthogonalize::Ortho4: Something's wrong. No splitting variables" );
	
		// select one randomly

		i = rand() % candidates2.size();
		prom = candidates2[i].var;
		ti1 = candidates2[i].ti;

		// split it

		if ( (*ti1).GetIM(prom) != '-' ) throw RuntimeException( "Orthogonalize4::Something's wrong - bad splitting variable" );
		newterm = *ti1;
		(*ti1).SetIM( prom, '1' );
		newterm.SetIM( prom, '0' );
		absorbed1 = RemoveAbsorbed( ti1 );
		absorbed2 = IsAbsorbed( newterm );
		if ( !absorbed2 ) ti2 = pla->AppendTerm( newterm );

			// recompute matrix
		tr1++;

		if ( !absorbed1 && absorbed2 ) {
			RecomputeMatrix();
			tr2++;
		}
		// else if ( !absorbed1 && !absorbed2 ) AppendMatrixRow( ti2 );   // this only slows it down
		else {
			GenerateMatrix();
			tr3++;
		}
		i = pla->GetTerms();
		MoveOKTerms();
		if ( pla->GetTerms() != i && pla->GetTerms() != 0 ) GenerateMatrix();
//		pla->SavePLA("bad.pla");
//		ok_terms->SavePLA("ok.pla");
	
	}

	if ( !IsOrthogonal(ok_terms) ) throw RuntimeException( "Orthogonalize4::Something's wrong - resulting PLA is not orthogonal" );
}

void Ortho::Orthogonalize5(void)
/*  1. find a term & variable to split, so that one (or both) split terms will be absorbed and, by this splitting,
		most intersected terms won't be intersected

*/
{
	Term newterm;
    list<Term>::iterator ti1, ti2, ti3, ti4;
    vector<int> candidates;
	vector<splitterm> candidates3;
	vector<splitterm> candidates2;
	vector<splitterm>::iterator cand3_it;
	vector<int> matrixrow;
    unsigned int prom;
	short pol;
    unsigned int i, j;
    list<Term> *data;
	unsigned int mx, cnt1, cnt0;
	unsigned int term;
	unsigned int val;
	bool absorbed1, absorbed2;		// are old/new terms absorbed?
	unsigned int orig_intersects, intersects;	// numbers of intersected terms
	bool consensus;				// consensus detected

	unsigned int tr1, tr2, tr3;
	tr1 = 0;
	tr2 = 0;
	tr3 = 0;

//	ofstream debugfile;
//	debugfile.open("ortho_debug.txt", std::ios::app);

	data = pla->GetData();
				// first preprocessing
	RemoveAbsorbed();
	pla->MergeTerms();
	MoveOKTerms2();
//	pla->SavePLA("bad.pla");
//	ok_terms->SavePLA("ok.pla");
	GenerateMatrix();

	while ( pla->GetTerms() != 0 ) {
			// find a term & variable to split, so that one (or both) split terms will be absorbed
		candidates3.clear();
		ti1 = data->begin();
		mx = 0;
		do {
			consensus = false;
			if ( (*ti1).GetDim() != 0 ) {
				orig_intersects = 0;		// compute # of intersecting terms
				for ( i = 0; i < pla->GetTerms(); i++ )
					if ( matrix[pla->GetIndex(ti1)][i] != -1 ) orig_intersects++;
				if ( orig_intersects >= mx ) {		// otherwise this term cannot compete
					for ( i = 0; i < pla->GetInputs(); i++ ) 
						if ( (*ti1).GetIM(i) == '-' ) {
							(*ti1).SetIM( i, '1' );
							absorbed1 = WillBeAbsorbed(ti1);
							(*ti1).SetIM( i, '0' );
							absorbed2 = WillBeAbsorbed(ti1);
							if ( absorbed1 && absorbed2 ) {		// it is the consenzus or could be merged. Remove it.
								// (*ti1).SetIM( i, '0' );		// is not necessary here. It was set above
								TryToMerge(ti1);
								(*ti1).SetIM( i, '1' );
								TryToMerge(ti1);
								cand3_it = candidates3.begin();
								while ( cand3_it != candidates3.end() ) {			// delete the term, if it was in candidates
									if ( (*cand3_it).ti == ti1 ) cand3_it = candidates3.erase(cand3_it);
									else cand3_it++;
								}
								ti1 = pla->DeleteTerm(ti1);
								GenerateMatrix();
								j = pla->GetTerms();
								MoveOKTerms();
								if ( pla->GetTerms() != j && pla->GetTerms() != 0 ) {
									GenerateMatrix();
									ti1 = data->begin();
									candidates3.clear();
								}
								//	ti1 = data->begin();		// this somehow produces worse results
								//	candidates3.clear();
								consensus = true;
								break;		// exit the for loop
							} else if ( absorbed1 ) {
								// (*ti1).SetIM( i, '0' );		// is not necessary here. It was set above
								intersects = 0;			// compute # of intersections
								for ( j = 0; j < pla->GetTerms(); j++ )
									if ( matrix[j][pla->GetIndex(ti1)] != -1 )
										if ( (*ti1).Intersects( *(pla->GetTermByIndex(j)) )) intersects++;
								if ( orig_intersects - intersects == mx && mx != 0 ) {			
									s_term.ti = ti1;
									s_term.var = i;
									s_term.pol = 0;
									candidates3.push_back(s_term);
								} else if ( orig_intersects - intersects > mx ) {			
									s_term.ti = ti1;
									s_term.var = i;
									s_term.pol = 0;
									candidates3.clear();
									candidates3.push_back(s_term);
									mx = orig_intersects - intersects;
								}
							} else if ( absorbed2 ) {
								(*ti1).SetIM( i, '1' );
								intersects = 0;			// compute # of intersections
								for ( j = 0; j < pla->GetTerms(); j++ )
									if ( matrix[j][pla->GetIndex(ti1)] != -1 )
										if ( (*ti1).Intersects( *(pla->GetTermByIndex(j)) )) intersects++;
								if ( orig_intersects - intersects == mx  && mx != 0 ) {			
									s_term.ti = ti1;
									s_term.var = i;
									s_term.pol = 1;
									candidates3.push_back(s_term);
								} else if ( orig_intersects - intersects > mx ) {			
									s_term.ti = ti1;
									s_term.var = i;
									s_term.pol = 1;
									candidates3.clear();
									candidates3.push_back(s_term);
									mx = orig_intersects - intersects;
								}
							}
							(*ti1).SetIM( i, '-' );
						}
				}
			}
			if (!consensus) ti1++;
		} while ( pla->GetTerms() != 0 && ti1 != data->end() );	  
		if ( pla->GetTerms() == 0 ) break;	// consensus eliminated all terms
		if ( !candidates3.empty() ) {				// one is absorbed
					// select one randomly
			i = rand() % candidates3.size();
			prom = candidates3[i].var;
			ti1 = candidates3[i].ti;
			pol = candidates3[i].pol;
			(*ti1).SetIM( prom, (1-pol)+'0' );	// try to merge the absorbed one
			TryToMerge(ti1);
			(*ti1).SetIM( prom, pol+'0' );
			RecomputeMatrix();
		} else
		{									// no chance for absobption
				// find term intersecting with most other terms
			mx = 0;
			candidates.clear();
			for ( ti1 = data->begin(); ti1 != data->end(); ti1++ ) 
				if ( (*ti1).GetDim() != 0 ) {
					orig_intersects = 0;
					for ( i = 0; i < pla->GetTerms(); i++ )
						if ( matrix[pla->GetIndex(ti1)][i] != -1 ) orig_intersects++;
					if ( orig_intersects == mx  && mx != 0 ) {
						candidates.push_back(pla->GetIndex(ti1));
					} else if ( orig_intersects > mx ) {			
						candidates.clear();
						candidates.push_back(pla->GetIndex(ti1));
						mx = orig_intersects;
					}
				}
				// select one randomly
			if ( candidates.empty() ) throw RuntimeException( "Orthogonalize::Ortho5: Something's wrong. No intersecting term." );
			term = candidates[rand() % candidates.size()];
			ti1 = pla->GetTermByIndex(term);

				// select the variable present in most intersected terms
			mx = 0;
			candidates.clear();
			for ( i = 0; i < pla->GetInputs(); i++ ) 
				if ( (*ti1).GetIM(i) == '-' ) {
					cnt0 = 0;
					cnt1 = 0;
					for ( ti2 = data->begin(); ti2 != data->end(); ti2++ ) {
						if ( matrix[term][pla->GetIndex(ti2)] != -1 )
							if ( (*ti2).GetIM(i) == '1' ) cnt1++;
							else if ( (*ti2).GetIM(i) == '0' ) cnt0++;
					}
					val = max( cnt0, cnt1 );
					if ( mx < val ) {
						candidates.clear();
						candidates.push_back(i);
						mx = val;
					} else if ( mx == val )
						candidates.push_back(i);
				}
			if ( candidates.empty() ) throw RuntimeException( "Orthogonalize::Ortho5: Something's wrong. No splitting variables" );
				// select one randomly
			prom = candidates[rand() % candidates.size()];

			newterm = *ti1;
			(*ti1).SetIM( prom, '1' );
			newterm.SetIM( prom, '0' );
			pla->AppendTerm( newterm );
			GenerateMatrix();
		}
		i = pla->GetTerms();
		RemoveAbsorbed();
		if ( pla->GetTerms() != i && pla->GetTerms() != 0 ) {
			GenerateMatrix();
			i = pla->GetTerms();
		}
		MoveOKTerms();
		if ( pla->GetTerms() != i && pla->GetTerms() != 0 ) GenerateMatrix();

		pla->SavePLA("bad.pla");
		ok_terms->SavePLA("ok.pla");
	}

	if ( !IsOrthogonal(ok_terms) ) throw RuntimeException( "Orthogonalize5::Something's wrong - resulting PLA is not orthogonal" );
}

bool Ortho::IsOrthogonal( PLA *_pla )
/*
    Checks if the PLA is orthogonal
*/
{
    list<Term>::const_iterator ti1, ti2;

    for ( ti1 = _pla->GetData()->begin(); ti1 != _pla->GetData()->end(); ti1++ )
        for ( ti2 = ti1; ti2 != _pla->GetData()->end(); ti2++ )
            if ( ti2 != ti1 )
                if ( (*ti1).Intersects(*ti2) )
					return false;
    return true;
}

PLA *Ortho::Run( const unsigned int method )
{
	switch (method) {
		case 1: Orthogonalize();
			break;
		case 2: Orthogonalize2();
			break;
		case 3: Orthogonalize3();
			break;
		case 4: Orthogonalize4();
			break;
		case 5: Orthogonalize5();
			break;
		default: throw RuntimeException( "Ortho:: Unknown orthogonalization method" );
	}

	ok_terms->SetNonUnique();
	ok_terms->MergeNeighbors();
	ok_terms->MergeTerms();

	if ( !IsOrthogonal( ok_terms ) ) throw RuntimeException( "Ortho::Something's wrong - resulting PLA is not orthogonal" );
	return ok_terms;
}
