/*Implementation of PStep_SeparationAlgorithm.h

It contains the member functions of PStep_SeparationAlgorithm5

*/

#include "PStep_SeparationAlgorithm.h"

//PStep_SeparationAlgorithm5:
//>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>

//Separate
//------------------------------------------------------------------------------------

//Separate capacity cuts using the algorithm of Lysgaard.
//These cuts are stored in AlgorithmMembers.CapacityCuts, new cuts are appended to the list.
//If it fails, separate homogeneous multistar cuts using the algorithm of Lysgaard.
//These cuts are stored in AlgorithmMembers.HomogeneousMultistarCuts, new cuts are appended to the list.
//If it fails, separate framed capacity cuts using the algorithm of Lysgaard.
//These cuts are stored in AlgorithmMembers.FramedCapacityCuts, new cuts are appended to the list.
//If that fails, separate strengthened comb cuts using the algorithm of Lysgaard.
//These cuts are stored in AlgorithmMembers.StrengthenedCombCuts, new cuts are appended to the list.
//If that also fails, separate 2 edge hypotour cuts using the algorithm of Lysgaard.
//These cuts are stored in AlgorithmMembers.TwoEdgeHypotourCuts, new cuts are appended to the list.
void PStep_SeparationAlgorithm5::Separate(){

	//Prepare the selected edges
	SetCurrentEdges();

	//Check whether there are any positively selected edges
	if (NumberOfEdgesPlusOne > 1){

		//Obtain current number of capacity cuts
		int NumberOfCCBeforeSeparation = AlgorithmMembers.NumberOfCapacityCuts;

		//Separate capacity cuts
		SeparateCCLysgaard();

		//Check whether new capacity cuts were added
		if (NumberOfCCBeforeSeparation == AlgorithmMembers.NumberOfCapacityCuts){

			//Obtain current number of framed capacity cuts
			int NumberOfFCCBeforeSeparation = AlgorithmMembers.NumberOfFramedCapacityCuts;

			//Separate framed capacity cuts
			SeparateFCCLysgaard();

			//Check whether new capacity cuts were added
			if (NumberOfFCCBeforeSeparation == AlgorithmMembers.NumberOfFramedCapacityCuts){
				//No violated capacity cuts were identified, separate homogeneous multistar cuts
				
				//Obtain current number of homogeneous mutlistar cuts
				int NumberOfHMCBeforeSeparation = AlgorithmMembers.NumberOfHomogeneousMultistarCuts;

				//Separate homogeneous multistar cuts
				SeparateHMCLysgaard();

				//Check whether new homogeneous multistar cuts were added
				if (NumberOfHMCBeforeSeparation == AlgorithmMembers.NumberOfHomogeneousMultistarCuts){
					//No violated homogeneous multistar cuts were identified, separate strengthened comb cuts

					//Obtain current number of strentghened comb cuts
					int NumberOfSCombBeforeSeparation = AlgorithmMembers.NumberOfSCombCuts;

					//Separate strengthened comb cuts
					SeparateSCombLysgaard();

					//Check whether new strengthened comb cuts were added
					if (NumberOfSCombBeforeSeparation == AlgorithmMembers.NumberOfSCombCuts){
						//No violated strengthened comb cuts were identified, separate 2 edge hypotour cuts

						//Separate 2 edge hypotour cuts
						Separate2EHLysgaard();

					}

				}

			}

		}

	}

}

//Homogeneous multistar cut separation functions:
//------------------------------------------------------------------------------------

//Separate homogeneus multistar cuts using the algorithm of Lysgaard.
//These cuts are stored in AlgorithmMembers.HomogeneousMultistar, new cuts are appended to the list.
void PStep_SeparationAlgorithm5::SeparateHMCLysgaard(){

	//Obtain start time of the separation algorithm
	timePoint tic = HomeMadeclock::now();
	
	//Initialize
	//-------------------------------------------------

	//Initialize support variables
	short int NumberOfCustomers = AlgorithmMembers.Instance.NumberOfCustomers;
	short int NumberOfLocations = NumberOfCustomers + 1;
	short int TotalNumberOfEdges = NumberOfLocations*(NumberOfLocations - 1) / 2;
	double Tolerance = AlgorithmMembers.Tolerance;
	std::vector < std::vector < PStep_Arc >> &Arcs = AlgorithmMembers.Arcs;
	double MinimumViolation = AlgorithmMembers.AlgorithmSettings.MinimumHMCViolation;

	//Obtain the maximum number of framed capacity cuts added per iteration
	int MaxCuts = AlgorithmMembers.AlgorithmSettings.MaxSCombPerIteration;

	//Initalize double
	double MaxViolation;
	//Initialize CnstrMgrPointer to contain the generated cuts
	CnstrMgrPointer MyCutsCMP;
	CMGR_CreateCMgr(&MyCutsCMP, MaxCuts);

	//Use Lysgaard's algorithm:
	//-------------------------------------------------
	MSTARSEP_SeparateMultiStarCuts(NumberOfCustomers, LysgaardDemand, int(AlgorithmMembers.Instance.VehicleCapacity), NumberOfEdgesPlusOne - 1, EdgeTail, EdgeHead, EdgeX, AlgorithmMembers.LysgaardCuts, MaxCuts, &MaxViolation, MyCutsCMP);

	//Process new cuts
	//-------------------------------------------------

	//Initialize indicators, indicating whether an edge is in use.
	//EdgeInUse[(j-1)*j/2+i] is corresponds to edge (i,j) for i<j.
	bool* EdgeInUse = new bool[TotalNumberOfEdges];
	//Initialize coefficients
	double* Coefficient = new double[TotalNumberOfEdges];

	//Go over all new cuts
	for (int i = 0; i<MyCutsCMP->Size; i++){

		//Initialize homogeneous multistar cut
		AlgorithmMembers.HomogeneousMultistarCuts.emplace_back();
		//Obtain new cut
		PStep_EdgeAndCoefficientCut &NewHomogeneousMultistarCut = AlgorithmMembers.HomogeneousMultistarCuts.back();
		//Update counter
		++AlgorithmMembers.NumberOfHomogeneousMultistarCuts;

		//Obtain coefficients A and B
		short int A = MyCutsCMP->CPL[i]->A;
		short int B = MyCutsCMP->CPL[i]->B;

		//Reset EdgeInUse and 
		for (short int j = 0; j < TotalNumberOfEdges; ++j){
			EdgeInUse[j] = false;
		}
		
		//Initialize indicators, indicating whether a location is in the Nucleus
		bool* LocationInNucleus = new bool[NumberOfLocations];
		for (short int j = 0; j < NumberOfLocations; ++j){
			LocationInNucleus[j] = false;
		}

		//Obtain Nucleus
		for (short int j = 1; j <= MyCutsCMP->CPL[i]->IntListSize; ++j){
			//Obtain current location
			short int k = MyCutsCMP->CPL[i]->IntList[j];
			if (k == NumberOfLocations){
				//Recode the depot
				k = 0;
			}
			//Add location the Nucleus
			LocationInNucleus[k] = true;
		}

		//Create list of locations that are in the Nucleus and not in the Nucleus
		std::list<short int> Nucleus;
		std::list<short int> NotInNucleus;
		//Go over all locations
		for (short int j = 0; j < NumberOfLocations; ++j){
			//Check whether the current location is in the Nucleus
			if (LocationInNucleus[j]){
				//Add the location to Nucleus
				Nucleus.push_back(j);
			}
			else{
				//Add the location to the NotInNucleus list
				NotInNucleus.push_back(j);
			}
		}

		//Set edge in use and coefficients of the edges out of the Nucleus:
		//Go over all locations in the Nucleus
		for (auto k = Nucleus.begin(); k != Nucleus.end(); ++k){
			//Go over all locations not in the Nucleus
			for (auto j = NotInNucleus.begin(); j != NotInNucleus.end(); ++j){
				//Obtain the edge index of the corresponding edge
				short int Index = Arcs[*k][*j].Index;
				//Check whether this edge is already in use
				if (EdgeInUse[Index]){
					//The edge is already in use, update the coefficient
					Coefficient[Index] += B;
				}
				else{
					//The edge is not yet in use, put it to use
					EdgeInUse[Index] = true;
					Coefficient[Index] = B;
				}
			}
		}

		//Initialize list of connectors
		std::list<short int> Connectors;
		//Obtain the connectors
		for (short int j = 1; j <= MyCutsCMP->CPL[i]->CListSize; j++){
			//Obtain current location
			short int k = MyCutsCMP->CPL[i]->CList[j];
			if (k == NumberOfLocations){
				//Recode the depot
				k = 0;
			}
			//Add connector to the list
			Connectors.push_back(k);
		}

		//Initialize list of satellites
		std::list<short int> Satellites;
		//Obtain the connectors
		for (short int j = 1; j <= MyCutsCMP->CPL[i]->ExtListSize; j++){
			//Obtain current location
			short int k = MyCutsCMP->CPL[i]->ExtList[j];
			if (k == NumberOfLocations){
				//Recode the depot
				k = 0;
			}
			//Add sattelite to the list
			Satellites.push_back(k);
		}

		//Set edge in use and coefficients of the edges between connectors and satellites:
		//Go over all connectors
		for (auto k = Connectors.begin(); k != Connectors.end(); ++k){
			//Go over all Sattelites
			for (auto j = Satellites.begin(); j != Satellites.end(); ++j){
				//Obtain the edge index of the corresponding edge
				short int Index = Arcs[*k][*j].Index;
				//Check whether this edge is already in use
				if (EdgeInUse[Index]){
					//The edge is already in use, update the coefficient
					Coefficient[Index] -= A;
				}
				else{
					//The edge is not yet in use, put it to use
					EdgeInUse[Index] = true;
					Coefficient[Index] = -A;
				}
			}
		}

		//Add all edges and coefficients to the cut
		for (short int j = 0; j < TotalNumberOfEdges; ++j){
			//Check whether the edge is in use
			if (EdgeInUse[j]){
				//The edge is in use, add it
				NewHomogeneousMultistarCut.EdgeAndCoefficients.emplace_back();
				NewHomogeneousMultistarCut.EdgeAndCoefficients.back().EdgeIndex = j;
				NewHomogeneousMultistarCut.EdgeAndCoefficients.back().Coefficient = Coefficient[j];
			}
		}

		//Add bound to the strengthened comb cut
		NewHomogeneousMultistarCut.Bound = MyCutsCMP->CPL[i]->L;

		//Check whether the violation of the cut is large enough
		if (ComputeViolation(NewHomogeneousMultistarCut, false) >= MinimumViolation){
			//The violation is large enough, so add the cut:

			//Move the current cut to LysgaardCuts
			CMGR_MoveCnstr(MyCutsCMP, AlgorithmMembers.LysgaardCuts, i, 0);
		}
		else{
			//Remove cut again
			AlgorithmMembers.HomogeneousMultistarCuts.pop_back();
			--AlgorithmMembers.NumberOfHomogeneousMultistarCuts;
		}

		//clean
		delete[]LocationInNucleus;

	}

	//clean
	delete[]EdgeInUse;
	delete[]Coefficient;
	CMGR_FreeMemCMgr(&MyCutsCMP);

	//Terminate:
	//-------------------------------------------------
	//Obtain end time of the separation algorithm
	timePoint toc = HomeMadeclock::now();
	//Update framed capacity cut separation time
	AlgorithmMembers.CurrentCGStatistics.TotalSeparateFCCTime += HomeMadeclock::seconds(tic, toc);

}