/* (C) Copyright 2003 Jens Lysgaard. All rights reserved. */
/* OSI Certified Open Source Software */
/* This software is licensed under the Common Public License Version 1.0 */

#include <stdlib.h>
#include <stdio.h>
#include "memmod.h"
#include "basegrph.h"
#include "strngcmp.h"
#include "cutbase.h"
#include "compress.h"

void COMPRESS_CheckV1Set(ReachPtr SupportPtr,
                         int NoOfCustomers,
                         int *CompNr,
                         double **XMatrix,
                         double *Slack,
                         int *CompListSize,
                         int *CompList,
                         int CutNr,
                         ReachPtr V1CutsPtr)
{
  int j,k;
  double XSumInSet;
  char *InNodeSet;
  char *CompInSet;

  InNodeSet = MemGetCV(NoOfCustomers+1);
  CompInSet = MemGetCV(NoOfCustomers+1);

  for (j=1; j<=NoOfCustomers; j++) InNodeSet[j] = 0;
  for (j=1; j<=NoOfCustomers; j++) CompInSet[j] = 0;

  *CompListSize = 0;

  for (k=1; k<=V1CutsPtr->LP[CutNr].CFN; k++)
  {
    j = V1CutsPtr->LP[CutNr].FAL[k];
    InNodeSet[j] = 1;

    if (CompInSet[CompNr[j]] == 0)
    {
      CompList[++(*CompListSize)] = CompNr[j];
      CompInSet[CompNr[j]] = 1;
    }
  }

  CUTBASE_CompXSumInSet(SupportPtr,NoOfCustomers,InNodeSet,NULL,0,
                        XMatrix,&XSumInSet);

  *Slack = (1.0 * V1CutsPtr->LP[CutNr].CFN) - XSumInSet - 1.0;

  MemFree(InNodeSet);
  MemFree(CompInSet);
}

void COMPRESS_ShrinkGraph(ReachPtr SupportPtr,
                          int NoOfCustomers,
                          double **XMatrix,
                          double **SMatrix,
                          int NoOfV1Cuts,
                          ReachPtr V1CutsPtr,
                          ReachPtr SAdjRPtr, /* Shrunk support graph */
                          ReachPtr SuperNodesRPtr, /* Nodes in supernodes */
                          int *ShrunkGraphCustNodes) /* Excl. depot */
{
  char NewLinks,TolerantShrinking,ShrinkableSet;
  int i,j,k,Tail,Head;
  int NoOfComponents,NodeListSize;
  int CompListSize;
  double EdgeEps,TripleEps,XVal;
  double MaxEdge;
  double Slack;
  double MaxShrinkTolerance,UsedTolerance;

  char *CVWrk1;
  int *IVWrk1, *IVWrk2, *IVWrk3, *IVWrk4;

  char *Shrinkable;
  int *CompNr; /* Node j is in component nr. CompNr[j]. */
  int *NodeList;
  int *CompList;

  ReachPtr CmprsEdgesRPtr;
  ReachPtr CompsRPtr;

  EdgeEps = 0.999; /* Shrink any edge >= Eps. */
  TripleEps = 1.999;

  TolerantShrinking = 0;
  MaxShrinkTolerance = 0.01;
  UsedTolerance = 0.0;

  CVWrk1 = MemGetCV(NoOfCustomers+1);

  IVWrk1 = MemGetIV(NoOfCustomers+1);
  IVWrk2 = MemGetIV(NoOfCustomers+1);
  IVWrk3 = MemGetIV(NoOfCustomers+1);
  IVWrk4 = MemGetIV(NoOfCustomers+1);

  Shrinkable = MemGetCV(NoOfCustomers+1);

  CompNr = MemGetIV(NoOfCustomers+2);
  NodeList = MemGetIV(NoOfCustomers+1);
  CompList = MemGetIV(NoOfCustomers+1);

  ReachInitMem(&CmprsEdgesRPtr,NoOfCustomers+1);
  ReachInitMem(&CompsRPtr,NoOfCustomers+1);

  for (i=1; i<NoOfCustomers; i++)
  {
    for (k=1; k<=SupportPtr->LP[i].CFN; k++)
    {
      j = SupportPtr->LP[i].FAL[k];
      if ((j <= NoOfCustomers) && (j > i))
      if (XMatrix[i][j] >= EdgeEps)
      {
        ReachAddForwArc(CmprsEdgesRPtr,i,j);
        ReachAddForwArc(CmprsEdgesRPtr,j,i);
      }
    }
  }

  do
  {
    ReachClearForwLists(CompsRPtr);

    ComputeStrongComponents(CmprsEdgesRPtr,CompsRPtr,
                            &NoOfComponents,NoOfCustomers,
                            CVWrk1,
                            IVWrk1,IVWrk2,IVWrk3,IVWrk4);

    for (i=1; i<=NoOfComponents; i++)
    {
      for (k=1; k<=CompsRPtr->LP[i].CFN; k++)
      {
        j = CompsRPtr->LP[i].FAL[k];
        CompNr[j] = i;
      }
    }

    CompNr[NoOfCustomers+1] = NoOfComponents+1; /* Depot in last comp. */

    for (i=1; i<=NoOfComponents+1; i++)
    for (j=1; j<=NoOfComponents+1; j++)
    SMatrix[i][j] = 0.0;

    for (i=1; i<=NoOfCustomers; i++) /* i = NoOfCustomers is ok, */
    {                                /* j may be the depot. */
      for (k=1; k<=SupportPtr->LP[i].CFN; k++)
      {
        j = SupportPtr->LP[i].FAL[k];
        if (j > i)
        {
          XVal = XMatrix[i][j];
          SMatrix[CompNr[i]][CompNr[j]] += XVal;

          if (CompNr[i] != CompNr[j])
          SMatrix[CompNr[j]][CompNr[i]] += XVal;
        }
      }
    }

    for (i=1; i<=NoOfComponents; i++)
    {
      if (SMatrix[i][i] < (CompsRPtr->LP[i].CFN - 1 + 0.01))
      Shrinkable[i] = 1;
      else
      Shrinkable[i] = 0;
    }

    /* Check for new possible compressions. */

    NewLinks = 0;

    for (i=1; i<NoOfComponents; i++)
    {
      if (Shrinkable[i] == 0) continue;

      for (j=i+1; j<=NoOfComponents; j++)
      {
        if (Shrinkable[j] == 0) continue;

        if (SMatrix[i][j] >= EdgeEps)
        {
          Tail = CompsRPtr->LP[i].FAL[1];
          Head = CompsRPtr->LP[j].FAL[1];

          ReachAddForwArc(CmprsEdgesRPtr,Tail,Head);
          ReachAddForwArc(CmprsEdgesRPtr,Head,Tail);
          NewLinks = 1;
        }
      }
    }

    if (NewLinks == 0)
    {
      /* (check triplets) */
      for (i=1; i<NoOfComponents; i++)
      {
        if (Shrinkable[i] == 0) continue;

        for (j=i+1; j<NoOfComponents; j++)
        {
          if (Shrinkable[j] == 0) continue;
          if (SMatrix[i][j] <= 0.01) continue;

          for (k=j+1; k<=NoOfComponents; k++)
          {
            if (Shrinkable[k] == 0) continue;
            if (SMatrix[i][k] <= 0.01) continue;
            if (SMatrix[j][k] <= 0.01) continue;

            XVal = SMatrix[i][j] + SMatrix[i][k] + SMatrix[j][k];

            if (XVal >= TripleEps)
            {
              Tail = CompsRPtr->LP[i].FAL[1];
              Head = CompsRPtr->LP[j].FAL[1];

              ReachAddForwArc(CmprsEdgesRPtr,Tail,Head);
              ReachAddForwArc(CmprsEdgesRPtr,Head,Tail);

              Head = CompsRPtr->LP[k].FAL[1];

              ReachAddForwArc(CmprsEdgesRPtr,Tail,Head);
              ReachAddForwArc(CmprsEdgesRPtr,Head,Tail);

              NewLinks = 1;
              goto EndOfTriplets;
            }
          }
        }
      }

      EndOfTriplets:

      i = 0;

    }

    if ((NewLinks == 0) &&
        (TolerantShrinking) &&
        (UsedTolerance < MaxShrinkTolerance))
    {
      Tail = 0;
      Head = 0;
      MaxEdge = 0.0;
      for (i=1; i<NoOfComponents; i++)
      for (j=i+1; j<=NoOfComponents; j++)
      {
        if (SMatrix[i][j] > MaxEdge)
        {
          MaxEdge = SMatrix[i][j];
          Tail = i;
          Head = j;
        }
      }

      if (Tail > 0)
      if ((1.0 - MaxEdge) <= (MaxShrinkTolerance - UsedTolerance))
      {
        NewLinks = 1;

        i = Tail;
        j = Head;

        Tail = CompsRPtr->LP[i].FAL[1];
        Head = CompsRPtr->LP[j].FAL[1];

        ReachAddForwArc(CmprsEdgesRPtr,Tail,Head);
        ReachAddForwArc(CmprsEdgesRPtr,Head,Tail);
        UsedTolerance += (1.0 - MaxEdge);
      }
    }

    if (NewLinks == 0)
    {
      for (i=1; i<=NoOfV1Cuts; i++)
      {
        COMPRESS_CheckV1Set(SupportPtr,NoOfCustomers,CompNr,
                            XMatrix,
                            &Slack,
                            &CompListSize,
                            CompList,
                            i,
                            V1CutsPtr);

        ShrinkableSet = 1;
        for (j=1; j<=CompListSize; j++)
        {
          if (Shrinkable[CompList[j]] == 0)
          {
            ShrinkableSet = 0;
            break;
          }
        }

        if ((Slack <= 0.01) && (CompListSize > 1) && (ShrinkableSet))
        {
          Tail = CompsRPtr->LP[CompList[1]].FAL[1];

          for (j=2; j<=CompListSize; j++)
          {
            Head = CompsRPtr->LP[CompList[j]].FAL[1];

            ReachAddForwArc(CmprsEdgesRPtr,Tail,Head);
            ReachAddForwArc(CmprsEdgesRPtr,Head,Tail);
          }

          NewLinks = 1;
          goto EndOfSetShrinking;
        }
      }

      EndOfSetShrinking:

      i = 0;
    }

  } while (NewLinks);

  /* SuperNode number NoOfComponents+1 is the depot. */

  for (i=1; i<=NoOfComponents+1; i++)
  {
    NodeListSize = 0;
    for (j=1; j<=NoOfComponents+1; j++)
    {
      if (j != i)
      if (SMatrix[i][j] >= 0.0001)
      {
        NodeList[++NodeListSize] = j;
      }
    }

    ReachSetForwList(SAdjRPtr,NodeList,i,NodeListSize);
  }

  for (i=1; i<=NoOfComponents; i++)
  {
    ReachSetForwList(SuperNodesRPtr,
                     CompsRPtr->LP[i].FAL,
                     i,
                     CompsRPtr->LP[i].CFN);
  }

  /* Depot: */
  NodeList[1] = NoOfCustomers+1;
  ReachSetForwList(SuperNodesRPtr,NodeList,NoOfComponents+1,1);

  *ShrunkGraphCustNodes = NoOfComponents;

  MemFree(CVWrk1);

  MemFree(IVWrk1);
  MemFree(IVWrk2);
  MemFree(IVWrk3);
  MemFree(IVWrk4);

  MemFree(CompNr);
  MemFree(NodeList);

  MemFree(Shrinkable);
  MemFree(CompList);

  ReachFreeMem(&CmprsEdgesRPtr);
  ReachFreeMem(&CompsRPtr);
}

