/* See {BigCube.h} */
#include <BigCube.h>

/* Last edited on 2007-01-27 14:44:08 by stolfi */

#define BigCube_C_COPYRIGHT \
  "Copyright © 2000,2007 Universidade Estadual de Campinas (UNICAMP)"

#include <Tools.h>
#include <Triangulation.h>
#include <Octf.h>
#include <OctfRep.h>
#include <OctfData.h>
#include <Squared.h>
#include <Cube.h>

#include <r4.h>
#include <vec.h>
#include <affirm.h>

#include <assert.h>
#include <stdio.h>

Place_t GlueBigCube(Place_t a, Place_t b, uint n);

Place_t GlueBigCube(Place_t a, Place_t b, uint n)
  { demand(n >= 1, "order must be at least 1");
    int M = 20;
    Place_t ta[M*M];
    Place_t tb[M*M];
    ta[0*M + 0] = a; tb[0*M + 0] = b;

    int i, j;
    for (i = 1; i < n; i++)
      { ta[i*M + 0] = Clock(PrevF(NextE(NextE(ta[(i-1)*M+0]))));
        tb[i*M + 0] = Clock(NextF(NextE(NextE(tb[(i-1)*M+0]))));
      }

    for (i = 0; i < n; i++)
      { for (j = 1; j < n; j++)
          { ta[i*M + j] = Clock(PrevE(PrevF(NextE(ta[i*M + (j-1)]))));
            tb[i*M + j] = Clock(PrevE(NextF(NextE(tb[i*M + (j-1)]))));
            assert(ta[i*M + j]!=a); 
            assert(tb[i*M + j]!=b);
          }
      }

    for (i = 0; i < n; i++)
      { Meld(tb[i*M + 0], ta[i*M + 0]);

        /* Update the edge slots */
        SetRingEdgeInfo(ta[i*M + 0], PWedge(ta[i*M + 0])->edge);
        SetRingEdgeInfo(NextE(ta[i*M + 0]), PWedge(NextE(ta[i*M + 0]))->edge);
        SetRingEdgeInfo(NextE(NextE(ta[i*M + 0])), PWedge(NextE(NextE(ta[i*M + 0])))->edge);
        SetRingEdgeInfo(NextE(NextE(NextE(ta[i*M + 0]))),
                   PWedge(NextE(NextE(NextE(ta[i*M + 0]))))->edge);

        /* Update the node slots */
        SetOrgAll(ta[i*M + 0], OrgV(ta[i*M + 0]));
        SetOrgAll(NextE(ta[i*M + 0]), OrgV(NextE(ta[i*M + 0])));
        SetOrgAll(NextE(NextE(ta[i*M + 0])), OrgV(NextE(NextE(ta[i*M + 0]))));
        SetOrgAll(NextE(NextE(NextE(ta[i*M + 0]))),
                   OrgV(NextE(NextE(NextE(ta[i*M + 0])))));
        SetOrgAll(Clock(ta[i*M + 0]), OrgV(Clock(ta[i*M + 0])));
        SetOrgAll(Clock(NextE(ta[i*M + 0])), OrgV(Clock(NextE(ta[i*M + 0]))));
        SetOrgAll(Clock(NextE(NextE(ta[i*M + 0]))),
                   OrgV(Clock(NextE(NextE(ta[i*M + 0])))));
        SetOrgAll(Clock(NextE(NextE(NextE(ta[i*M + 0])))),
                   OrgV(Clock(NextE(NextE(NextE(ta[i*M + 0]))))));

        SetPneg(ta[i*M + 0], PnegP(tb[i*M + 0]));
        SetPneg(PrevE(ta[i*M + 0]), PnegP(PrevE(tb[i*M + 0])));
        SetPneg(PrevE(PrevE(ta[i*M + 0])), PnegP(PrevE(PrevE(tb[i*M + 0]))));
        SetPneg(PrevE(PrevE(PrevE(ta[i*M + 0]))),
                PnegP(PrevE(PrevE(PrevE(tb[i*M + 0])))));

      }
    
    for (i = 0; i < n; i++)
      { for (j = 1; j < n; j++)
          { Meld(tb[i*M + j],ta[i*M + j]);

            /* Update the edge slots */
            SetRingEdgeInfo(ta[i*M + j], PWedge(ta[i*M + j])->edge);
            SetRingEdgeInfo(NextE(ta[i*M + j]), PWedge(NextE(ta[i*M + j]))->edge);
            SetRingEdgeInfo(NextE(NextE(ta[i*M + j])),
                                 PWedge(NextE(NextE(ta[i*M + j])))->edge);
            SetRingEdgeInfo(NextE(NextE(NextE(ta[i*M + j]))),
                       PWedge(NextE(NextE(NextE(ta[i*M + j]))))->edge);

            /* Update the node slots */
            SetOrgAll(ta[i*M + j], OrgV(ta[i*M + j]));
            SetOrgAll(NextE(ta[i*M + j]), OrgV(NextE(ta[i*M + j])));
            SetOrgAll(NextE(NextE(ta[i*M + j])), OrgV(NextE(NextE(ta[i*M + j]))));
            SetOrgAll(NextE(NextE(NextE(ta[i*M + j]))),
                                       OrgV(NextE(NextE(NextE(ta[i*M + j])))));
            SetOrgAll(Clock(ta[i*M + j]), OrgV(Clock(ta[i*M + j])));
            SetOrgAll(Clock(NextE(ta[i*M + j])), OrgV(Clock(NextE(ta[i*M + j]))));
            SetOrgAll(Clock(NextE(NextE(ta[i*M + j]))),
                                         OrgV(Clock(NextE(NextE(ta[i*M + j])))));
            SetOrgAll(Clock(NextE(NextE(NextE(ta[i*M + j])))),
                       OrgV(Clock(NextE(NextE(NextE(ta[i*M + j]))))));

            SetPneg(ta[i*M + j], PnegP(tb[i*M + j]));
            SetPneg(PrevE(ta[i*M + j]), PnegP(PrevE(tb[i*M + j])));
            SetPneg(PrevE(PrevE(ta[i*M + j])), PnegP(PrevE(PrevE(tb[i*M + j]))));
            SetPneg(PrevE(PrevE(PrevE(ta[i*M + j]))),
                    PnegP(PrevE(PrevE(PrevE(tb[i*M + j])))));                 
          }
      }

    return ta[0*M + 0];
  }

FourPlaces_vec_t Make1DCubearray(uint order)
  { SixPlaces_vec_t ca = SixPlaces_vec_new(order);
    FourPlaces_vec_t cb = FourPlaces_vec_new(order);
    int i, j;
    for (i = 0; i < order; i++)
      { ca.el[i] = MakeCube();
        cb.el[i].p[0] = ca.el[i].p[1]; 
        cb.el[i].p[1] = ca.el[i].p[3]; 
        cb.el[i].p[2] = ca.el[i].p[0]; 
        cb.el[i].p[3] = ca.el[i].p[5];
      }

    /* Glue the cubes: */
    for (j = 0; j < order-1; j++)
      { GlueCube(ca.el[j].p[2], PrevF(ca.el[j+1].p[3])); }
    
    return cb;
  }

TwoPlaces_vec_t Make2DCubeArray(uint order)
  { FourPlaces_vec_t cb = FourPlaces_vec_new(order* order);
    TwoPlaces_vec_t cc = TwoPlaces_vec_new(order*order);

    /* Assembling the array {cb} */
    int i,j,k;
    for (i = 0; i < order; i++)
      { FourPlaces_vec_t ca = Make1DCubearray(order);
        for (j = 0; j < order; j++)
          { for (k = 0; k < 4; k++)
              { cb.el[i*order + j].p[k] = ca.el[j].p[k]; }
          }
      }

    /* Gluing the rows: */
    for (j = 0; j < (order-1); j++)
      { for (k = 0; k < order; k++)
          { GlueCube(cb.el[j*order + k].p[1], Clock(NextE(NextE(cb.el[(j+1)*order + k].p[0])))); }
      }

    /* Selecting the places to return: */
    for (i = 0; i < order; i++)
      { for (j = 0; j < order; j++)
        { for (k = 0; k < 2; k++)
            { cc.el[i*order + j].p[k] = cb.el[i*order + j].p[k+2]; }
        }
      }
      
    return cc;
  }

TwoPlaces_vec_t MakeBigCube(uint order)
  { FourPlaces_vec_t cd = FourPlaces_vec_new(order*order*order);
    TwoPlaces_vec_t ce = TwoPlaces_vec_new(order*order);
    int i,j,k,l;
    for (i = 0; i < order; i++)
      { TwoPlaces_vec_t cc = Make2DCubeArray(order);
        for (j = 0; j < order; j++)
          { for (k = 0; k < order; k++)
              { for (l = 0; l < 2; l++)
                  { cd.el[(i*order + j)*order + k].p[l] = cc.el[j*order + k].p[l]; }
              }
          }
      }

    /* Gluing: */
    for (i = 0; i < (order-1); i++)
      { for (k = 0; k < order; k++)
          { for (j = 0; j < order; j++)
              { Place_t pikj = cd.el[(i*order + k)*order +j].p[1];
                Place_t qikj = cd.el[((i+1)*order +k)*order + j].p[0];
                GlueCube(pikj, Clock(qikj));
              }
          }
      }

    /* Selecting the places for return */
    for (j = 0; j < order; j++)
      { for (k = 0; k < order; k++)
          { ce.el[j*order + k].p[0] = cd.el[(0*order + k)*order + j].p[0];
            ce.el[j*order + k].p[1] = cd.el[((order-1)*order + k)*order + j].p[1];
          }
      }
    return ce;
  }

void SetCubeProperties(TwoPlaces_vec_t *a, uint n, Topology_t *tp)
  { uint index = 0;
    Place_vec_t c1 = Place_vec_new(n);
    Place_vec_t c2 = Place_vec_new(n);
    Place_vec_t c3 = Place_vec_new(n);
    Place_vec_t c4 = Place_vec_new(n);

    /* Set the original nodes */

    SetExNode(Clock(NextE(a->el[0*n + 0].p[0])), FALSE);    /* node 1 */   
    SetExNode(Clock(a->el[(n-1)*n + 0].p[0]), FALSE);         /* node 2 */
    SetExNode(a->el[(n-1)*n + (n-1)].p[0], FALSE);              /* node 3 */
    SetExNode(a->el[(n-1)*n + 0].p[1], FALSE);                /* node 4 */
    SetExNode(Clock(a->el[(n-1)*n + (n-1)].p[1]), FALSE);       /* node 5 */
    SetExNode(Clock(NextE(a->el[0*n + (n-1)].p[1])), FALSE);  /* node 6 */
    SetExNode(PrevE(a->el[0*n + 0].p[1]), FALSE);           /* node 7 */
    SetExNode(PrevE(a->el[0*n + (n-1)].p[0]), FALSE);         /* node 8 */

    /* Set the original edgess */
    int i, j, k;
    for (i = 0; i < n; i++)
      { SetExEdge(PrevE(a->el[i*n + 0].p[1]), FALSE);
        SetExEdge(NextE(a->el[i*n + 0].p[0]), FALSE);
        SetExEdge(a->el[(n-1)*n + i].p[0], FALSE);
        SetExEdge(a->el[(n-1)*n + i].p[1], FALSE);
        SetExEdge(PrevE(a->el[(i)*n + (n-1)].p[0]), FALSE);
        SetExEdge(NextE(a->el[(i)*n + (n-1)].p[1]), FALSE);
        SetExEdge(NextE(NextE(a->el[0*n + i].p[1])), FALSE);
        SetExEdge(NextE(NextE(a->el[0*n + i].p[0])), FALSE);
      }

    /* computing @{edge->?}s on the colums */
    c1.el[0] = NextE(NextF(a->el[(n-1)*n + 0].p[0]));
    c2.el[0] = PrevE(NextF(a->el[(n-1)*n + (n-1)].p[0]));
    c3.el[0] = NextE(NextF(NextE(a->el[0*n + 0].p[0])));
    c4.el[0] = PrevE(NextF(PrevE(a->el[0*n + (n-1)].p[0])));

    for (k = 1; k < n; k++)
      { c1.el[k] = Clock(PrevE(PrevF(NextE(c1.el[k-1]))));
        c2.el[k] = Clock(NextE(PrevF(PrevE(c2.el[k-1]))));
        c3.el[k] = Clock(PrevE(PrevF(NextE(c3.el[k-1]))));
        c4.el[k] = Clock(NextE(PrevF(PrevE(c4.el[k-1]))));
      }

    /* Now set the original edgess */
    for (j = 0; j < n; j++)
      { SetExEdge(c1.el[j], FALSE);
        SetExEdge(c2.el[j], FALSE);
        SetExEdge(c3.el[j], FALSE);
        SetExEdge(c4.el[j], FALSE);
      }

    for (i = 0; i < tp->NF; i++)
      { Place_t f = tp->wall.el[i]->pa;
          { if (WallOnBoundary(f)){ SetExWall(f, FALSE); index++; } }
      }
    assert(index == 6*n*n);
  } 

TwoPlaces_vec_t BuildBigRawCube(uint order, bool_t original)
  { TwoPlaces_vec_t cd = TwoPlaces_vec_new(order*order*order);
    int i,j,k,l;
    for (i = 0; i < order; i++)
      { TwoPlaces_vec_t cc = Make2DCubeArray(order);
        for (j = 0; j < order; j++)
          { for (k = 0; k < order; k++)
              { for (l = 0; l < 2; l++)
                  { cd.el[(i*order + j)*order + k].p[l] = cc.el[j*order + k].p[l]; }
              }
          }
      }

    /* Gluing: */
    for (i = 0; i < (order-1); i++)
      { for (k = 0; k < order; k++)
          { for (j = 0; j < order; j++)
              { Place_t pikj = cd.el[(i*order + k)*order + j].p[1];
                Place_t qikj = cd.el[((i+1)*order + k)*order + j].p[0];
                GlueCube(pikj, Clock(qikj));
              }
          }
      }
    
    return cd;
  }

#define BigCube_C_author \
  "Created by L. A. P. Lozada, 2000.\n" \
  "Modification history:\n" \
  "  30-08-2000 : [???] Nice version of the {MakeOctahedron} procedure.\n" \
  "  19-09-2000 : [???] Added the procedure "MakeDodecahedronTriang".\n" \
  "  27-10-2000 : [???] Modified {SetCubeProperties} procedure.\n" \
  "  04-11-2000 : [???] Added the {CubeNegNodes} and {CubeBarycenter} procedures.\n" \
  "  26-01-2007 : [???] Moved {MakeBigRawCube} from {MakeBigCubeFixed.c} to here."