/*Implementation of PStep_SeparationAlgorithm.h

It contains the member functions of PStep_SeparationAlgorithm1

*/

#include "PStep_SeparationAlgorithm.h"

//PStep_SeparationAlgorithm1:
//>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>

//Constructor functions
//------------------------------------------------------------------------------------

//Basic constructor
PStep_SeparationAlgorithm1::PStep_SeparationAlgorithm1(PStep_AlgorithmMembers& InitAlgorithmMembers)
	:AlgorithmMembers(InitAlgorithmMembers), EdgeHead(NULL), EdgeTail(NULL), EdgeX(NULL)
{
	//Obtain number of locations
	short int NumberOfLocations = AlgorithmMembers.Instance.NumberOfCustomers +1;
	//Obtain Demand
	double* Demand = AlgorithmMembers.Instance.Demand;

	//Obtain demand as required by Lysgaard
	//Lysgaard requires a 0 as first entry, and Demand[i] to correspond to customer i.
	//Also demand is required to be integer.
	LysgaardDemand = new int[NumberOfLocations];
	LysgaardDemand[0] = 0;
	for (short int i = 1; i<NumberOfLocations; i++){
		LysgaardDemand[i] = int(Demand[i-1]);
	}

}

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

//Separate capacity cuts using the algorithm of Lysgaard.
//The cuts are stored in AlgorithmMembers.CapacityCuts, new cuts are appended to the list.
void PStep_SeparationAlgorithm1::Separate(){

	//Prepare the selected edges
	SetCurrentEdges();

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

		//Separate capacity cuts
		SeparateCCLysgaard();

	}

}

//Capacity constraint separation functions:
//------------------------------------------------------------------------------------

//Separate capacity constraints using the package of Lysgaard.
//The cuts are stored in AlgorithmMembers.CapacityCuts, new cuts are appended to the list.
void PStep_SeparationAlgorithm1::SeparateCCLysgaard(){

	//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;
	double Tolerance = AlgorithmMembers.Tolerance;
	double MinimumViolation = AlgorithmMembers.AlgorithmSettings.MinimumCCViolation;

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

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

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

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

		//Initialize list that temporarily stores all customers in the current capacity cut.
		//The value 1 corresponds to the first customer.
		std::list<short int> Customers;

		//Obtain all customers in the subset.
		for (int j = 1; j <= MyCutsCMP->CPL[i]->IntListSize; j++){
			//Add customer to subset
			short int CurrentCustomer = MyCutsCMP->CPL[i]->IntList[j];
			Customers.push_back(CurrentCustomer);
		}

		//Convert the list of customers and total demand to an EdgeFlowCut:
		//Create capacity cut
		PStep_EdgeFlowCut NewCapacityCut;

		//Go over all locations in the cut
		for (auto CurrentLocationIndex = Customers.begin(); CurrentLocationIndex != Customers.end(); ++CurrentLocationIndex){
			//Initialize support variables:
			//Obtain the current locations
			short int CurrentLocation = *CurrentLocationIndex;
			//Obtain the first next index
			std::list<short int>::iterator FirstNextIndex = CurrentLocationIndex;
			++FirstNextIndex;
			//Go over all other locations in the cut
			for (auto OtherLocationIndex = FirstNextIndex; OtherLocationIndex != Customers.end(); ++OtherLocationIndex){
				//Obtain the other location
				short int OtherLocation = *OtherLocationIndex;
				//Add the edge between the current and other location to the cut
				if (CurrentLocation < OtherLocation){
					NewCapacityCut.Indices.push_back((OtherLocation - 1)*OtherLocation / 2 + CurrentLocation);
				}
				else{
					NewCapacityCut.Indices.push_back((CurrentLocation - 1)*CurrentLocation / 2 + OtherLocation);
				}
			}
		}
		//Add bound to the capacity cut
		NewCapacityCut.Bound = MyCutsCMP->CPL[i]->RHS;

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

			//Add the new capacity constraint to the set of capacity constraints
			AlgorithmMembers.CapacityCuts.push_back(NewCapacityCut);
			//Update capacity cut counter
			++AlgorithmMembers.NumberOfCapacityCuts;
			//Move the current cut to LysgaardCuts
			CMGR_MoveCnstr(MyCutsCMP, AlgorithmMembers.LysgaardCuts, i, 0);

		}

	}

	//clean
	CMGR_FreeMemCMgr(&MyCutsCMP);

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

}

//Other functions:
//------------------------------------------------------------------------------------

//Set the separation members pertaining to the used edges in the current solution,
//i.e., NumberOfEdgesPlusOne, EdgeHead, EdgeTail, EdgeX.
void PStep_SeparationAlgorithm1::SetCurrentEdges(){

	//Reset edge information
	NumberOfEdgesPlusOne = 1;
	delete[]EdgeHead;
	EdgeHead = NULL;
	delete[]EdgeTail;
	EdgeTail = NULL;
	delete[]EdgeX;
	EdgeX = NULL;

	//Obtain support variables
	short int NumberOfLocations = AlgorithmMembers.Instance.NumberOfCustomers + 1;
	double Tolerance = AlgorithmMembers.Tolerance;

	//Obtain current edge fraction
	double* EdgeFraction = AlgorithmMembers.CurrentSolution.EdgeFraction;

	//Identify positively selected edges:
	//-------------------------------------------------

	//Initialize vectors containing the positively selected edges
	std::vector<int> EdgeTailVector;
	std::vector<int> EdgeHeadVector;
	std::vector<double> EdgeXVector;
	//Go over all edges
	int Counter = 0;
	for (int j = 1; j < NumberOfLocations; ++j){
		//check whether the edge fraction from the depot (i=0) is positive
		if (EdgeFraction[Counter] > Tolerance){
			//Update number of positive edges
			++NumberOfEdgesPlusOne;
			//Store edge information
			EdgeTailVector.push_back(j);
			EdgeHeadVector.push_back(NumberOfLocations); //NumberOfCustomers+1 is used by Lysgaard to represent the depot.
			EdgeXVector.push_back(EdgeFraction[Counter]);
		}
		//Update Counter
		++Counter;
		//Go over the other locations
		for (int i = 1; i < j; ++i){
			//check whether the edge fraction is positive
			if (EdgeFraction[Counter] > Tolerance){
				//Update number of positive edges
				++NumberOfEdgesPlusOne;
				//Store edge information
				EdgeTailVector.push_back(i);
				EdgeHeadVector.push_back(j);
				EdgeXVector.push_back(EdgeFraction[Counter]);
			}
			//Update counter
			++Counter;
		}
	}

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

		//Construct input for the algorithm by Lysgaard:
		//-------------------------------------------------

		//Initialize arrays with edge information.
		EdgeTail = new int[NumberOfEdgesPlusOne];
		EdgeHead = new int[NumberOfEdgesPlusOne];
		EdgeX = new double[NumberOfEdgesPlusOne];
		//Lysgaard requires a bogus first edge in the array.
		EdgeTail[0] = 0;
		EdgeHead[0] = 0;
		EdgeX[0] = 0;
		//Go over all identified edges and add them to the arrays
		for (int i = 1; i < NumberOfEdgesPlusOne; ++i){
			EdgeTail[i] = EdgeTailVector[i - 1];
			EdgeHead[i] = EdgeHeadVector[i - 1];
			EdgeX[i] = EdgeXVector[i - 1];
		}
	}

}

//Compute the amount by which a constraint is violated.
//The output is non-positive if there is no violation.
//The values of the edge selection are taken from AlgorithmMembers.CurrentSolution.EdgeFraction.
//If bound of CurrentEdgeFlowCut acts as upper bound, set UpperBound to true, otherwise false.
double PStep_SeparationAlgorithm1::ComputeViolation(PStep_EdgeFlowCut &CurrentEdgeFlowCut, bool UpperBound){

	//Obtain EdgeFraction
	double* (&EdgeFraction) = AlgorithmMembers.CurrentSolution.EdgeFraction;

	//Initialize violation as minus the bound
	double Violation = CurrentEdgeFlowCut.Bound;

	//Go over all edges
	for (auto CurrentEdge = CurrentEdgeFlowCut.Indices.begin(); CurrentEdge != CurrentEdgeFlowCut.Indices.end(); ++CurrentEdge){
		Violation -= EdgeFraction[*CurrentEdge];
	}

	//Output violation
	if (UpperBound){
		return -Violation;
	}
	return Violation;

}

//Compute the amount by which a constraint is violated.
//The output is non-positive if there is no violation.
//The values of the edge selection are taken from AlgorithmMembers.CurrentSolution.EdgeFraction.
//If bound of PStep_EdgeAndCoefficientCut acts as upper bound, set UpperBound to true, otherwise false.
double PStep_SeparationAlgorithm1::ComputeViolation(PStep_EdgeAndCoefficientCut &CurrentEdgeAndCoefficientCut, bool UpperBound){

	//Obtain EdgeFraction
	double* (&EdgeFraction) = AlgorithmMembers.CurrentSolution.EdgeFraction;

	//Initialize violation as minus the bound
	double Violation = CurrentEdgeAndCoefficientCut.Bound;

	//Go over all edges
	for (auto CurrentEdgeAndCoefficient = CurrentEdgeAndCoefficientCut.EdgeAndCoefficients.begin(); CurrentEdgeAndCoefficient != CurrentEdgeAndCoefficientCut.EdgeAndCoefficients.end(); ++CurrentEdgeAndCoefficient){
		Violation -= (CurrentEdgeAndCoefficient->Coefficient)*(EdgeFraction[CurrentEdgeAndCoefficient->EdgeIndex]);
	}

	//Output violation
	if (UpperBound){
		return -Violation;
	}
	return Violation;

}

//Destructor functions
//------------------------------------------------------------------------------------

//Default destructor
PStep_SeparationAlgorithm1::~PStep_SeparationAlgorithm1(){

	delete[]LysgaardDemand;
	delete[]EdgeHead;
	delete[]EdgeTail;
	delete[]EdgeX;

}