/*************************************************************************
 * Copyright (c) 2011 AT&T Intellectual Property
 * All rights reserved. This program and the accompanying materials
 * are made available under the terms of the Eclipse Public License v1.0
 * which accompanies this distribution, and is available at
 * https://www.eclipse.org/legal/epl-v10.html
 *
 * Contributors: Details at https://graphviz.org
 *************************************************************************/

/* tlayout.c:
 * Written by Emden R. Gansner
 *
 * Module for initial layout, using point nodes and ports.
 *
 * Note: If interior nodes are not connected, they tend to fly apart,
 * despite being tied to port nodes. This occurs because, as initially
 * coded, as the nodes tend to straighten into a line, the radius
 * grows causing more expansion. Is the problem really here and not
 * with disconnected nodes in xlayout? If here, we can either forbid
 * expansion or eliminate repulsion between nodes only connected
 * via port nodes.
 */

#include "config.h"

/* uses PRIVATE interface */
#define FDP_PRIVATE 1

#include <math.h>
#include <stdbool.h>
#include <stdlib.h>
#include <sys/types.h>
#include <time.h>
#ifndef _WIN32
#include <unistd.h>
#endif
#include <fdpgen/dbg.h>
#include <fdpgen/grid.h>
#include <neatogen/neato.h>

#ifndef HAVE_SRAND48
#define srand48 srand
#endif

#include <common/globals.h>
#include <fdpgen/tlayout.h>

#define D_useGrid (fdp_parms->useGrid)
#define D_useNew (fdp_parms->useNew)
#define D_numIters (fdp_parms->numIters)
#define D_unscaled (fdp_parms->unscaled)
#define D_C (fdp_parms->C)
#define D_Tfact (fdp_parms->Tfact)
#define D_K (fdp_parms->K)
#define D_T0 (fdp_parms->T0)

/* Actual parameters used; initialized using fdp_parms, then possibly
 * updated with graph-specific values.
 */
typedef struct {
  int useGrid;  /* use grid for speed up */
  int useNew;   /* encode x-K into attractive force */
  long seed;    /* seed for position RNG */
  int numIters; /* actual iterations in layout */
  int maxIters; /* max iterations in layout */
  int unscaled; /* % of iterations used in pass 1 */
  double C;     /* Repulsion factor in xLayout */
  double Tfact; /* scale temp from default expression */
  double K;     /* spring constant; ideal distance */
  double T0;    /* initial temperature */
  int smode;    /* seed mode */
  double Cell;  /* grid cell size */
  double Wd;    /* half-width of boundary */
  double Ht;    /* half-height of boundary */
  int pass1;    /* iterations used in pass 1 */
  int loopcnt;  /* actual iterations in this pass */
} parms_t;

static parms_t parms;

#define T_useGrid (parms.useGrid)
#define T_useNew (parms.useNew)
#define T_seed (parms.seed)
#define T_numIters (parms.numIters)
#define T_maxIters (parms.maxIters)
#define T_unscaled (parms.unscaled)
#define T_C (parms.C)
#define T_Tfact (parms.Tfact)
#define T_K (parms.K)
#define T_T0 (parms.T0)
#define T_smode (parms.smode)
#define T_Cell (parms.Cell)
#define T_Wd (parms.Wd)
#define T_Ht (parms.Ht)
#define T_pass1 (parms.pass1)
#define T_loopcnt (parms.loopcnt)

#define EXPFACTOR 1.2
#define DFLT_maxIters 600
#define DFLT_K 0.3
#define DFLT_Cell 0.0
#define DFLT_seed 1
#define DFLT_smode INIT_RANDOM

static double cool(int t) { return T_T0 * (T_maxIters - t) / T_maxIters; }

static void reset_params(void) { T_T0 = -1.0; }

/* Set parameters for expansion phase based on initial
 * layout parameters. If T0 is not set, we set it here
 * based on the size of the graph. In this case, we
 * return true, so that fdp_tLayout can unset T0, to be
 * reset by a recursive call to fdp_tLayout.
 *
 * @return Does `reset` need to be called later?
 */
static bool init_params(graph_t *g, xparams *xpms) {
  bool ret = false;

  if (T_T0 == -1.0) {
    int nnodes = agnnodes(g);

    T_T0 = T_Tfact * T_K * sqrt(nnodes) / 5;
#ifdef DEBUG
    if (Verbose) {
      prIndent();
      fprintf(stderr, "tlayout %s", agnameof(g));
      fprintf(stderr, "(%s) : T0 %f\n", agnameof(GORIG(g->root)), T_T0);
    }
#endif
    ret = true;
  }

  xpms->T0 = cool(T_pass1);
  xpms->K = T_K;
  xpms->C = T_C;
  xpms->numIters = T_maxIters - T_pass1;

  if (T_numIters >= 0) {
    if (T_numIters <= T_pass1) {
      T_loopcnt = T_numIters;
      xpms->loopcnt = 0;
    } else if (T_numIters <= T_maxIters) {
      T_loopcnt = T_pass1;
      xpms->loopcnt = T_numIters - T_pass1;
    }
  } else {
    T_loopcnt = T_pass1;
    xpms->loopcnt = xpms->numIters;
  }
  return ret;
}

/// initialize parameters based on root graph attributes
void fdp_initParams(graph_t *g) {
  T_useGrid = D_useGrid;
  T_useNew = D_useNew;
  T_numIters = D_numIters;
  T_unscaled = D_unscaled;
  T_Cell = DFLT_Cell;
  T_C = D_C;
  T_Tfact = D_Tfact;
  T_maxIters =
      late_int(g, agattr_text(g, AGRAPH, "maxiter", NULL), DFLT_maxIters, 0);
  D_K = T_K = late_double(g, agattr_text(g, AGRAPH, "K", NULL), DFLT_K, 0.0);
  if (D_T0 == -1.0) {
    T_T0 = late_double(g, agattr_text(g, AGRAPH, "T0", NULL), -1.0, 0.0);
  } else
    T_T0 = D_T0;
  T_seed = DFLT_seed;
  T_smode = setSeed(g, DFLT_smode, &T_seed);
  if (T_smode == INIT_SELF) {
    agwarningf("fdp does not support start=self - ignoring\n");
    T_seed = DFLT_smode;
  }

  T_pass1 = T_unscaled * T_maxIters / 100;

  if (T_useGrid) {
    if (T_Cell <= 0.0)
      T_Cell = 3 * T_K;
  }
#ifdef DEBUG
  if (Verbose) {
    prIndent();
    fprintf(stderr, "Params %s : K %f T0 %f Tfact %f maxIters %d unscaled %d\n",
            agnameof(g), T_K, T_T0, T_Tfact, T_maxIters, T_unscaled);
  }
#endif
}

static void doRep(node_t *p, node_t *q, double xdelta, double ydelta,
                  double dist2) {
  double force;
  double dist;

  while (dist2 == 0.0) {
    xdelta = 5 - rand() % 10;
    ydelta = 5 - rand() % 10;
    dist2 = xdelta * xdelta + ydelta * ydelta;
  }
  if (T_useNew) {
    dist = sqrt(dist2);
    force = T_K * T_K / (dist * dist2);
  } else
    force = T_K * T_K / dist2;
  if (IS_PORT(p) && IS_PORT(q))
    force *= 10.0;
  DISP(q)[0] += xdelta * force;
  DISP(q)[1] += ydelta * force;
  DISP(p)[0] -= xdelta * force;
  DISP(p)[1] -= ydelta * force;
}

/// repulsive force = K × K ÷ d or K × K ÷ d × d
static void applyRep(Agnode_t *p, Agnode_t *q) {
  double xdelta, ydelta;

  xdelta = ND_pos(q)[0] - ND_pos(p)[0];
  ydelta = ND_pos(q)[1] - ND_pos(p)[1];
  doRep(p, q, xdelta, ydelta, xdelta * xdelta + ydelta * ydelta);
}

static void doNeighbor(Grid *grid, int i, int j, node_list *nodes) {
  cell *cellp = findGrid(grid, i, j);
  node_list *qs;
  Agnode_t *p;
  Agnode_t *q;
  double xdelta, ydelta;
  double dist2;

  if (cellp) {
#ifdef DEBUG
    if (Verbose >= 3) {
      prIndent();
      fprintf(stderr, "  doNeighbor (%d,%d) : %d\n", i, j, gLength(cellp));
    }
#endif
    for (; nodes != 0; nodes = nodes->next) {
      p = nodes->node;
      for (qs = cellp->nodes; qs != 0; qs = qs->next) {
        q = qs->node;
        xdelta = (ND_pos(q))[0] - (ND_pos(p))[0];
        ydelta = (ND_pos(q))[1] - (ND_pos(p))[1];
        dist2 = xdelta * xdelta + ydelta * ydelta;
        if (dist2 < T_Cell * T_Cell)
          doRep(p, q, xdelta, ydelta, dist2);
      }
    }
  }
}

static int gridRepulse(void *c, void *g) {
  cell *const cellp = c;
  Grid *const grid = g;
  node_list *nodes = cellp->nodes;
  int i = cellp->p.i;
  int j = cellp->p.j;
  node_list *p;
  node_list *q;

#ifdef DEBUG
  if (Verbose >= 3) {
    prIndent();
    fprintf(stderr, "gridRepulse (%d,%d) : %d\n", i, j, gLength(cellp));
  }
#endif
  for (p = nodes; p != 0; p = p->next) {
    for (q = nodes; q != 0; q = q->next)
      if (p != q)
        applyRep(p->node, q->node);
  }

  doNeighbor(grid, i - 1, j - 1, nodes);
  doNeighbor(grid, i - 1, j, nodes);
  doNeighbor(grid, i - 1, j + 1, nodes);
  doNeighbor(grid, i, j - 1, nodes);
  doNeighbor(grid, i, j + 1, nodes);
  doNeighbor(grid, i + 1, j - 1, nodes);
  doNeighbor(grid, i + 1, j, nodes);
  doNeighbor(grid, i + 1, j + 1, nodes);

  return 0;
}

/* Attractive force = weight × (d × d) ÷ K
 *  or        force = (d - L(e)) × weight(e)
 */
static void applyAttr(Agnode_t *p, Agnode_t *q, Agedge_t *e) {
  double xdelta, ydelta;
  double force;
  double dist;
  double dist2;

  xdelta = ND_pos(q)[0] - ND_pos(p)[0];
  ydelta = ND_pos(q)[1] - ND_pos(p)[1];
  dist2 = xdelta * xdelta + ydelta * ydelta;
  while (dist2 == 0.0) {
    xdelta = 5 - rand() % 10;
    ydelta = 5 - rand() % 10;
    dist2 = xdelta * xdelta + ydelta * ydelta;
  }
  dist = sqrt(dist2);
  if (T_useNew)
    force = ED_factor(e) * (dist - ED_dist(e)) / dist;
  else
    force = ED_factor(e) * dist / ED_dist(e);
  DISP(q)[0] -= xdelta * force;
  DISP(q)[1] -= ydelta * force;
  DISP(p)[0] += xdelta * force;
  DISP(p)[1] += ydelta * force;
}

static void updatePos(Agraph_t *g, double temp, bport_t *pp) {
  Agnode_t *n;
  double temp2;
  double len2;
  double x, y, d;
  double dx, dy;

  temp2 = temp * temp;
  for (n = agfstnode(g); n; n = agnxtnode(g, n)) {
    if (ND_pinned(n) & P_FIX)
      continue;
    dx = DISP(n)[0];
    dy = DISP(n)[1];
    len2 = dx * dx + dy * dy;

    /* limit by temperature */
    if (len2 < temp2) {
      x = ND_pos(n)[0] + dx;
      y = ND_pos(n)[1] + dy;
    } else {
      double fact = temp / sqrt(len2);
      x = ND_pos(n)[0] + dx * fact;
      y = ND_pos(n)[1] + dy * fact;
    }

    /* if ports, limit by boundary */
    if (pp) {
      d = sqrt(x * x / (T_Wd * T_Wd) + y * y / (T_Ht * T_Ht));
      if (IS_PORT(n)) {
        ND_pos(n)[0] = x / d;
        ND_pos(n)[1] = y / d;
      } else if (d >= 1.0) {
        ND_pos(n)[0] = 0.95 * x / d;
        ND_pos(n)[1] = 0.95 * y / d;
      } else {
        ND_pos(n)[0] = x;
        ND_pos(n)[1] = y;
      }
    } else {
      ND_pos(n)[0] = x;
      ND_pos(n)[1] = y;
    }
  }
}

#define FLOOR(d) ((int)floor(d))

static void gAdjust(Agraph_t *g, double temp, bport_t *pp, Grid *grid) {
  Agnode_t *n;
  Agedge_t *e;

  if (temp <= 0.0)
    return;

  clearGrid(grid);

  for (n = agfstnode(g); n; n = agnxtnode(g, n)) {
    DISP(n)[0] = DISP(n)[1] = 0;
    addGrid(grid, FLOOR((ND_pos(n))[0] / T_Cell),
            FLOOR((ND_pos(n))[1] / T_Cell), n);
  }

  for (n = agfstnode(g); n; n = agnxtnode(g, n)) {
    for (e = agfstout(g, n); e; e = agnxtout(g, e))
      if (n != aghead(e))
        applyAttr(n, aghead(e), e);
  }
  walkGrid(grid, gridRepulse);

  updatePos(g, temp, pp);
}

static void adjust(Agraph_t *g, double temp, bport_t *pp) {
  Agnode_t *n;
  Agnode_t *n1;
  Agedge_t *e;

  if (temp <= 0.0)
    return;

  for (n = agfstnode(g); n; n = agnxtnode(g, n)) {
    DISP(n)[0] = DISP(n)[1] = 0;
  }

  for (n = agfstnode(g); n; n = agnxtnode(g, n)) {
    for (n1 = agnxtnode(g, n); n1; n1 = agnxtnode(g, n1)) {
      applyRep(n, n1);
    }
    for (e = agfstout(g, n); e; e = agnxtout(g, e)) {
      if (n != aghead(e))
        applyAttr(n, aghead(e), e);
    }
  }

  updatePos(g, temp, pp);
}

/* Create initial layout of nodes
 * TODO :
 *  Position nodes near neighbors with positions.
 *  Use bbox to reset K.
 */
static pointf initPositions(graph_t *g, bport_t *pp) {
  int nG = agnnodes(g) - NPORTS(g);
  double size;
  Agnode_t *np;
  int n_pos = 0; /* no. of nodes with position info */
  boxf bb = {{0, 0}, {0, 0}};
  pointf ctr; /* center of boundary ellipse */
  long local_seed;
  double PItimes2 = M_PI * 2.0;

  for (np = agfstnode(g); np; np = agnxtnode(g, np)) {
    if (ND_pinned(np)) {
      if (n_pos) {
        bb.LL.x = MIN(ND_pos(np)[0], bb.LL.x);
        bb.LL.y = MIN(ND_pos(np)[1], bb.LL.y);
        bb.UR.x = MAX(ND_pos(np)[0], bb.UR.x);
        bb.UR.y = MAX(ND_pos(np)[1], bb.UR.y);
      } else {
        bb.UR.x = bb.LL.x = ND_pos(np)[0];
        bb.UR.y = bb.LL.y = ND_pos(np)[1];
      }
      n_pos++;
    }
  }

  size = T_K * (sqrt((double)nG) + 1.0);
  T_Wd = T_Ht = EXPFACTOR * (size / 2.0);
  if (n_pos == 1) {
    ctr.x = bb.LL.x;
    ctr.y = bb.LL.y;
  } else if (n_pos > 1) {
    double alpha, area, width, height, quot;
    ctr.x = (bb.LL.x + bb.UR.x) / 2.0;
    ctr.y = (bb.LL.y + bb.UR.y) / 2.0;
    width = EXPFACTOR * (bb.UR.x - bb.LL.x);
    height = EXPFACTOR * (bb.UR.y - bb.LL.y);
    area = 4.0 * T_Wd * T_Ht;
    quot = width * height / area;
    if (quot >= 1.0) { /* If bbox has large enough area, use it */
      T_Wd = width / 2.0;
      T_Ht = height / 2.0;
    } else if (quot > 0.0) { /* else scale up to have enough area */
      quot = 2.0 * sqrt(quot);
      T_Wd = width / quot;
      T_Ht = height / quot;
    } else { /* either width or height is 0 */
      if (width > 0) {
        height = area / width;
        T_Wd = width / 2.0;
        T_Ht = height / 2.0;
      } else if (height > 0) {
        width = area / height;
        T_Wd = width / 2.0;
        T_Ht = height / 2.0;
      }
      /* If width = height = 0, use Wd and Ht as defined above for
       * the case the n_pos == 0.
       */
    }

    /* Construct enclosing ellipse */
    alpha = atan2(T_Ht, T_Wd);
    T_Wd = T_Wd / cos(alpha);
    T_Ht = T_Ht / sin(alpha);
  } else {
    ctr.x = ctr.y = 0;
  }

  /* Set seed value */
  if (T_smode == INIT_RANDOM)
    local_seed = T_seed;
  else {
#if defined(_WIN32)
    local_seed = (long)time(NULL);
#else
    local_seed = getpid() ^ time(NULL);
#endif
  }
  srand48(local_seed);

  /* If ports, place ports on and nodes within an ellipse centered at origin
   * with halfwidth Wd and halfheight Ht.
   * If no ports, place nodes within a rectangle centered at origin
   * with halfwidth Wd and halfheight Ht. Nodes with a given position
   * are translated. Wd and Ht are set to contain all positioned points.
   * The reverse translation will be applied to all
   * nodes at the end of the layout.
   * TODO: place unfixed points using adjacent ports or fixed pts.
   */
  if (pp) {
    while (pp->e) { /* position ports on ellipse */
      np = pp->n;
      ND_pos(np)[0] = T_Wd * cos(pp->alpha) + ctr.x;
      ND_pos(np)[1] = T_Ht * sin(pp->alpha) + ctr.y;
      ND_pinned(np) = P_SET;
      pp++;
    }
    for (np = agfstnode(g); np; np = agnxtnode(g, np)) {
      if (IS_PORT(np))
        continue;
      if (ND_pinned(np)) {
        ND_pos(np)[0] -= ctr.x;
        ND_pos(np)[1] -= ctr.y;
      } else {
        pointf p = {0.0, 0.0};
        int cnt = 0;
        node_t *op;
        edge_t *ep;
        for (ep = agfstedge(g, np); ep; ep = agnxtedge(g, ep, np)) {
          if (aghead(ep) == agtail(ep))
            continue;
          op = aghead(ep) == np ? agtail(ep) : aghead(ep);
          if (!hasPos(op))
            continue;
          if (cnt) {
            p.x = (p.x * cnt + ND_pos(op)[0]) / (cnt + 1);
            p.y = (p.y * cnt + ND_pos(op)[1]) / (cnt + 1);
          } else {
            p.x = ND_pos(op)[0];
            p.y = ND_pos(op)[1];
          }
          cnt++;
        }
        if (cnt > 1) {
          ND_pos(np)[0] = p.x;
          ND_pos(np)[1] = p.y;
        } else if (cnt == 1) {
          ND_pos(np)[0] = 0.98 * p.x + 0.1 * ctr.x;
          ND_pos(np)[1] = 0.9 * p.y + 0.1 * ctr.y;
        } else {
          double angle = PItimes2 * drand48();
          double radius = 0.9 * drand48();
          ND_pos(np)[0] = radius * T_Wd * cos(angle);
          ND_pos(np)[1] = radius * T_Ht * sin(angle);
        }
        ND_pinned(np) = P_SET;
      }
    }
  } else {
    if (n_pos) { /* If positioned nodes */
      for (np = agfstnode(g); np; np = agnxtnode(g, np)) {
        if (ND_pinned(np)) {
          ND_pos(np)[0] -= ctr.x;
          ND_pos(np)[1] -= ctr.y;
        } else {
          ND_pos(np)[0] = T_Wd * (2.0 * drand48() - 1.0);
          ND_pos(np)[1] = T_Ht * (2.0 * drand48() - 1.0);
        }
      }
    } else { /* No ports or positions; place randomly */
      for (np = agfstnode(g); np; np = agnxtnode(g, np)) {
        ND_pos(np)[0] = T_Wd * (2.0 * drand48() - 1.0);
        ND_pos(np)[1] = T_Ht * (2.0 * drand48() - 1.0);
      }
    }
  }

  return ctr;
}

/* Given graph g with ports nodes, layout g respecting ports.
 * If some node have position information, it may be useful to
 * reset temperature and other parameters to reflect this.
 */
void fdp_tLayout(graph_t *g, xparams *xpms) {
  bport_t *const pp = PORTS(g);
  const bool reset = init_params(g, xpms);
  const pointf ctr = initPositions(g, pp);

  if (T_useGrid) {
    Grid *const grid = mkGrid(agnnodes(g));
    adjustGrid(grid, agnnodes(g));
    for (int i = 0; i < T_loopcnt; i++) {
      const double temp = cool(i);
      gAdjust(g, temp, pp, grid);
    }
    delGrid(grid);
  } else {
    for (int i = 0; i < T_loopcnt; i++) {
      const double temp = cool(i);
      adjust(g, temp, pp);
    }
  }

  if (ctr.x != 0.0 || ctr.y != 0.0) {
    for (Agnode_t *n = agfstnode(g); n; n = agnxtnode(g, n)) {
      ND_pos(n)[0] += ctr.x;
      ND_pos(n)[1] += ctr.y;
    }
  }
  if (reset)
    reset_params();
}