#include "GR_assert.h"
#include "GR_Geometry.h"
#include "GR_Mesh2D.h"
#include "GR_Vec.h"
#include "GR_events.h"
#include "GR_BFace.h"

//@ Insert a point on a boundary face
int Mesh2D::
iInsertOnBdryFace (Vert * const pVNew, Face * const pF,
		   Cell * const pC, const bool qSwap, const bool qForce)
{
  int iRetVal = 0;
  if (!qAllowBdryChanges)
    return -1;

  // Logging state on exit should be INSERT_VERTEX if it starts out that way. 
  bool qExitAsInsertVertex = false;
#ifdef SUMAA_LOG
  if (SUMAA_LOG_current == INSERT_VERTEX) {
    qExitAsInsertVertex = true;
  }
  else {
    qExitAsInsertVertex = false;
    SUMAA_LOG_EVENT_BEGIN(INSERT_VERTEX);
  }
  assert (SUMAA_LOG_current == INSERT_VERTEX);
#endif
  assert (pF->iFullCheck ());
  assert (pC->iFullCheck ());
  assert (pF->pCCellLeft () == pC || pF->pCCellRight () == pC);
  assert (pC->eType () == Cell::eTriCell);
  assert (pF->pCCellOpposite (pC)->eType () == Cell::eBdryEdge);

  TriCell *pTC = dynamic_cast<TriCell *>(pC);
  BFace *pBF = dynamic_cast<BFace *>(pF->pCCellOpposite (pC));

  Vert *pVTwin;
  BdryPatch *pBP = pBF->pPatchPointer();
  if (pBP->qIsPeriodic()) {
    // Identify the twin bdry face and split it.
    Periodic *pP = pBP->pPGetPeriodicData();
    Face *pFTwin = const_cast<Face*>(pP->pFFaceOpposite(pBP, pF));
    assert(pFTwin);
    BFace *pBFTwin = dynamic_cast<BFace*>(pFTwin->pCCellLeft());
    if (!pBFTwin) pBFTwin = dynamic_cast<BFace*>(pFTwin->pCCellRight());
    assert(pBFTwin); // Better be able to find the twin!

    pVTwin = pVNewVert();
    (*pVTwin) = *pVNew;
    if (qForce) {
      // Map coordinates
      double adTwinCoords[2];
      pP->vPeriodicLocation(pVNew->adCoords(), pBP, adTwinCoords);
      pVTwin->vSetCoords(2, adTwinCoords);
    }

    Cell *pCTwin = pF->pCCellOpposite(pBFTwin);
    iRetVal += iInsertOnBdryFace(pVTwin, pFTwin, pCTwin, qSwap, qForce);

    // If there's some way to at least register verts as being on a
    // periodic bdry, then it'll be easier to efficiently find twin
    // faces / edges.  Not entirely sure yet where this data should be
    // stored, though...
    // vRegisterPeriodicVert(pVNew, pP, pVTwin);
    // vRegisterPeriodicVert(pVTwin, pP, pVNew);
  }
  //  Won't need this since the info is kept in a patch now...
  //  int iBC = pBF->iBdryCond();

  // Should now calculate the location of the new vertex.
  // Use the EdgeBFace->vGetBisectionPoint for that...

  //  fprintf(stdout, "Before call to vGetBisectionPoint()\n");
  double adMid[2];
  if (!qForce) {
    dynamic_cast<BdryEdge *>(pBF)->vGetBisectionPoint (adMid);

#ifndef NDEBUG
//      fprintf (stdout, "The new coords are (%f, %f)\n", adMid[0], adMid[1]);
#endif
    pVNew->vSetCoords (2, adMid);
  }

  Vert *pVA, *pVB, *pVC;
  if (pF->pCCellLeft () == pC) {
    pVA = pF->pVVert (0);
    pVB = pF->pVVert (1);
  }
  else {
    pVB = pF->pVVert (0);
    pVA = pF->pVVert (1);
  }
  pVC = pTC->pVVertOpposite (pF);

  Face *pFA = pTC->pFFaceOpposite (pVA);
  Face *pFB = pTC->pFFaceOpposite (pVB);
  Face *pFC = pTC->pFFaceOpposite (pVC);

  assert (iOrient2D (pVA, pVNew, pVC) == 1);
  assert (iOrient2D (pVC, pVNew, pVB) == 1);

  Cell *pC1 = pC;
  Cell *pC2 = pCNewCell ();
  // Make sure that both cells are in the same region as the original
  pC1->vSetRegion (pC->iRegion ());
  pC2->vSetRegion (pC->iRegion ());

  BFace *pBF1 = pBF;
  BFace *pBF2 = pBFNewBFace ();

  Face *pFC1 = pFC;
  Face *pFC2 = pFNewFace ();
  Face *pFD = pFNewFace ();

  //@@ Fix up new edges
  pFC1->vAssign (pC1, pBF1, pVNew, pVB);
  pFC1->vSetFaceLoc (Face::eBdryFace);
  pFC2->vAssign (pC2, pBF2, pVA, pVNew);
  pFC2->vSetFaceLoc (Face::eBdryFace);
  pFD->vAssign (pC1, pC2, pVC, pVNew);

  //@@ Adjust one old edge
  pFB->vReplaceCell (pC, pC2);

  pBF2->vAssign (pFC2);
  // Make sure patch is set too...
  dynamic_cast<BdryEdge *>(pBF2)->
    vSetPatch (dynamic_cast<BdryEdge *>(pBF)->iPatchIndex (),
	       dynamic_cast<BdryEdge *>(pBF)->pBdryPointer ());

  //@@ Fix up cells
  pC1->vAssign (pFA, pFC1, pFD);
  pC2->vAssign (pFB, pFC2, pFD);

  //@@ Add a vertex hint for the new vert
  pVNew->vSetType (Vert::eBdry);
  pVNew->vSetHintFace (pFD);
  pVA->vSetHintFace (pFC2);
  pVB->vSetHintFace (pFC1);

  assert (pC1->iFullCheck ());
  assert (pC2->iFullCheck ());
  assert (pFA->iFullCheck ());
  assert (pFB->iFullCheck ());
  assert (pFC1->iFullCheck ());
  assert (pFC2->iFullCheck ());
  assert (pFD->iFullCheck ());

  //@@ Swap locally
  SUMAA_LOG_EVENT_END (INSERT_VERTEX);
  SUMAA_LOG_EVENT_BEGIN (NET_SWAPPING);
  if (qSwap) {
    // If you're going to swap while inserting, all recurse and use
    // the Delaunay criterion.
    iRetVal += iFaceSwap (pFA);
    iRetVal += iFaceSwap (pFB);
  }
  SUMAA_LOG_EVENT_END (NET_SWAPPING);
  if (qExitAsInsertVertex) {
    SUMAA_LOG_EVENT_BEGIN (INSERT_VERTEX);
  }
  spVUpdateLS.insert(pVNew);
  return (iRetVal);
}

//@ Insert a point on a face separating two regions in the mesh.
int Mesh2D::
iInsertOnInternalBdryFace (Vert * const pVNew, Face * const pF,
			   const bool qSwap, const bool qForce)
  // Note: Geometrically, pVNew may not fall exactly on pF because
  // of round-off.  Once things are done, though, the connectivity will
  // be the same as if it had.
{
  //  Before:               After:
  //              B                    B
  //             /|\                  /|\         .
  //         A1 / | \ A2          A1 / | \ A2
  //           /  |  \              /  MB \       .
  //          /   |   \            /1  |  2\      .
  //         /    |    \          /_N1_N_N2_\     .
  //        C\  1 | 2  /D        C\    |    /D
  //          \   |   /            \3  |  4/
  //           \  |  /              \  MA /
  //         B1 \ | / B2          B1 \ | / B2
  //             \|/                  \|/
  //              A                    A
  //

  if (!qAllowBdryChanges)
    return -1;

  double adMid[2];

  // Logging state on exit should be INSERT_VERTEX if it starts out that way.
  bool qExitAsInsertVertex;
#ifdef SUMAA_LOG
  if (SUMAA_LOG_current == INSERT_VERTEX) {
    qExitAsInsertVertex = true;
  }
  else {
    qExitAsInsertVertex = false;
    SUMAA_LOG_EVENT_BEGIN(INSERT_VERTEX);
  }
  assert (SUMAA_LOG_current == INSERT_VERTEX);
#endif
  assert (pF->iFullCheck ());
  assert (pF->iFaceLoc () == Face::eBdryTwoSide);

  // The coords for pVNew should be calculated using the information
  // stored in the patches...

  // Will pF really be a pointer to a boundary face...? Let's hope so.
  if (!qForce) {
    dynamic_cast<BdryEdge *>(pF)->vGetBisectionPoint (adMid);

#ifndef NDEBUG
//      fprintf (stdout, "New coords: (%f, %f)\n",
//	     adMid[0], adMid[1]);
#endif

    // Set the coords
    pVNew->vSetCoords (2, adMid);
  }

  Cell *pC1 = pF->pCCellLeft ();
  Cell *pC2 = pF->pCCellRight ();

  assert (pC1->iFullCheck ());
  assert (pC2->iFullCheck ());
  assert (pC1->iRegion () != pC2->iRegion ());
  assert (pC1->eType () == Cell::eTriCell);
  assert (pC2->eType () == Cell::eTriCell);

  TriCell *const pTC1 = dynamic_cast<TriCell *>(pC1);
  TriCell *const pTC2 = dynamic_cast<TriCell *>(pC2);

  Vert *pVA = pF->pVVert (0);
  Vert *pVB = pF->pVVert (1);
  Vert *pVC = pTC1->pVVertOpposite (pF);
  Vert *pVD = pTC2->pVVertOpposite (pF);

  assert (iOrient2D (pVA, pVB, pVC) == 1);
  assert (iOrient2D (pVB, pVA, pVD) == 1);

  Face *pFA1 = pTC1->pFFaceOpposite (pVA);
  Face *pFB1 = pTC1->pFFaceOpposite (pVB);
  Face *pFA2 = pTC2->pFFaceOpposite (pVA);
  Face *pFB2 = pTC2->pFFaceOpposite (pVB);
  Face *pFM = pF;

  assert (iOrient2D (pVA, pVD, pVNew) == 1);
  assert (iOrient2D (pVD, pVB, pVNew) == 1);
  assert (iOrient2D (pVB, pVC, pVNew) == 1);
  assert (iOrient2D (pVC, pVA, pVNew) == 1);

  //@@ Create new cells
  Cell *pC3 = pCNewCell ();
  Cell *pC4 = pCNewCell ();
  pC3->vSetRegion (pC1->iRegion ());
  pC4->vSetRegion (pC2->iRegion ());

  //@@ Create and set up new faces
  Face *pFMA = pFNewFace ();
  Face *pFMB = pFM;
  Face *pFN1 = pFNewFace ();
  Face *pFN2 = pFNewFace ();
  pFMA->vSetFaceLoc (Face::eBdryTwoSide);
  pFMB->vSetFaceLoc (Face::eBdryTwoSide);
  pFN1->vSetFaceLoc (Face::eInterior);
  pFN2->vSetFaceLoc (Face::eInterior);

  pFMA->vAssign (pC3, pC4, pVA, pVNew);
  pFMB->vAssign (pC1, pC2, pVNew, pVB);
  pFN1->vAssign (pC1, pC3, pVC, pVNew);
  pFN2->vAssign (pC4, pC2, pVD, pVNew);

  //@@ Fix two outer faces which have new cells
  pFB1->vReplaceCell (pC1, pC3);
  pFB2->vReplaceCell (pC2, pC4);

  //@@ Set up the new cells
  pC1->vAssign (pFA1, pFN1, pFMB);
  pC2->vAssign (pFA2, pFN2, pFMB);
  pC3->vAssign (pFB1, pFN1, pFMA);
  pC4->vAssign (pFB2, pFN2, pFMA);

  //@@ Add a vertex hint for the new vert
  pVNew->vSetType (Vert::eBdryTwoSide);
  pVNew->vSetHintFace (pFN1);
  pVA->vSetHintFace (pFMA);
  pVB->vSetHintFace (pFMB);

  //@@ Swap locally
  SUMAA_LOG_EVENT_END (INSERT_VERTEX);
  SUMAA_LOG_EVENT_BEGIN (NET_SWAPPING);
  int iRetVal = 0;
  if (qSwap) {
    iRetVal += iFaceSwap (pFA1);
    iRetVal += iFaceSwap (pFA2);
    iRetVal += iFaceSwap (pFB1);
    iRetVal += iFaceSwap (pFB2);
  }
  SUMAA_LOG_EVENT_END (NET_SWAPPING);
  if (qExitAsInsertVertex) {
    SUMAA_LOG_EVENT_BEGIN (INSERT_VERTEX);
  }
  spVUpdateLS.insert(pVNew);
  return (iRetVal);
}

//@ Insert a point on an interior face
int Mesh2D::
iInsertOnFace (Vert * const pVNew, Face * const pF,
	       Cell * const pC, const bool qSwap, const bool qForce)
  // Note: Geometrically, pVNew may not fall exactly on pF because
  // of round-off.  Once things are done, though, the connectivity will
  // be the same as if it had.
{
  assert (pF->iFullCheck ());
  Cell *pC1 = pF->pCCellLeft ();
  Cell *pC2 = pF->pCCellRight ();
  assert (pC1->iFullCheck ());
  assert (pC2->iFullCheck ());
  assert (pC1 == pC || pC2 == pC);

  if (pC1->eType () == Cell::eBdryEdge) {
    iOnBdryEdge2++;
    return iInsertOnBdryFace (pVNew, pF, pC2, qSwap, qForce);
  }
  if (pC2->eType () == Cell::eBdryEdge) {
    iOnBdryEdge2++;
    return iInsertOnBdryFace (pVNew, pF, pC1, qSwap, qForce);
  }
  if (pF->iFaceLoc () == Face::eBdryTwoSide) {
    iInternalBdry2++;
    return iInsertOnInternalBdryFace (pVNew, pF, qSwap, qForce);
  }
  assert (pC1->iRegion () == pC2->iRegion ());

  // Logging state on exit should be INSERT_VERTEX if it starts out that way.
  bool qExitAsInsertVertex;
#ifdef SUMAA_LOG
  if (SUMAA_LOG_current == INSERT_VERTEX) {
    qExitAsInsertVertex = true;
  }
  else {
    qExitAsInsertVertex = false;
    SUMAA_LOG_EVENT_BEGIN(INSERT_VERTEX);
  }
  assert (SUMAA_LOG_current == INSERT_VERTEX);
#endif
  assert (pC->iFullCheck ());

  iOnEdge++;

  assert (pC1->eType () == Cell::eTriCell);
  assert (pC2->eType () == Cell::eTriCell);

  TriCell *const pTC1 = dynamic_cast<TriCell *>(pC1);
  TriCell *const pTC2 = dynamic_cast<TriCell *>(pC2);

  Vert *pVA = pF->pVVert (0);
  Vert *pVB = pF->pVVert (1);
  Vert *pVC = pTC1->pVVertOpposite (pF);
  Vert *pVD = pTC2->pVVertOpposite (pF);

  assert (iOrient2D (pVA, pVB, pVC) == 1);
  assert (iOrient2D (pVB, pVA, pVD) == 1);

  Face *pFA1 = pTC1->pFFaceOpposite (pVA);
  Face *pFB1 = pTC1->pFFaceOpposite (pVB);
  Face *pFA2 = pTC2->pFFaceOpposite (pVA);
  Face *pFB2 = pTC2->pFFaceOpposite (pVB);
  Face *pFM = pF;

  assert (iOrient2D (pVA, pVD, pVNew) == 1);
  assert (iOrient2D (pVD, pVB, pVNew) == 1);
  assert (iOrient2D (pVB, pVC, pVNew) == 1);
  assert (iOrient2D (pVC, pVA, pVNew) == 1);

  //@@ Create new cells
  Cell *pC3 = pCNewCell ();
  Cell *pC4 = pCNewCell ();
  pC3->vSetRegion (pC1->iRegion ());
  pC4->vSetRegion (pC1->iRegion ());

  //@@ Create and set up new faces
  Face *pFMA = pFNewFace ();
  Face *pFMB = pFM;
  Face *pFN1 = pFNewFace ();
  Face *pFN2 = pFNewFace ();
  pFMA->vSetFaceLoc (Face::eInterior);
  pFN1->vSetFaceLoc (Face::eInterior);
  pFN2->vSetFaceLoc (Face::eInterior);

  pFMA->vAssign (pC3, pC4, pVA, pVNew);
  pFMB->vAssign (pC1, pC2, pVNew, pVB);
  pFN1->vAssign (pC1, pC3, pVC, pVNew);
  pFN2->vAssign (pC4, pC2, pVD, pVNew);

  //@@ Fix two outer faces which have new cells
  pFB1->vReplaceCell (pC1, pC3);
  pFB2->vReplaceCell (pC2, pC4);

  //@@ Set up the new cells
  pC1->vAssign (pFA1, pFN1, pFMB);
  pC2->vAssign (pFA2, pFN2, pFMB);
  pC3->vAssign (pFB1, pFN1, pFMA);
  pC4->vAssign (pFB2, pFN2, pFMA);

  //@@ Add a vertex hint for the new vert
  pVNew->vSetType (Vert::eInterior);
  pVNew->vSetHintFace (pFN1);
  pVA->vSetHintFace (pFMA);
  pVB->vSetHintFace (pFMB);

  //@@ Swap locally
  SUMAA_LOG_EVENT_END (INSERT_VERTEX);
  SUMAA_LOG_EVENT_BEGIN (NET_SWAPPING);
  int iRetVal = 0;
  if (qSwap) {
    iRetVal += iFaceSwap (pFA1);
    iRetVal += iFaceSwap (pFA2);
    iRetVal += iFaceSwap (pFB1);
    iRetVal += iFaceSwap (pFB2);
  }
  SUMAA_LOG_EVENT_END (NET_SWAPPING);
  if (qExitAsInsertVertex) {
    SUMAA_LOG_EVENT_BEGIN (INSERT_VERTEX);
  }
  spVUpdateLS.insert(pVNew);
  return (iRetVal);
}

//@ Insert a point inside a cell
int Mesh2D::
iInsertInInterior (Vert * const pVNew, Cell * const pC,
		   const bool qSwap)
{
  // Logging state on exit should be INSERT_VERTEX if it starts out that way.
  bool qExitAsInsertVertex;
#ifdef SUMAA_LOG
  if (SUMAA_LOG_current == INSERT_VERTEX) {
    qExitAsInsertVertex = true;
  }
  else {
    qExitAsInsertVertex = false;
    SUMAA_LOG_EVENT_BEGIN(INSERT_VERTEX);
  }
  assert (SUMAA_LOG_current == INSERT_VERTEX);
#endif
  assert (pC->iFullCheck ());
  assert (pC->eType () == Cell::eTriCell);
  TriCell *pTC = dynamic_cast<TriCell *>(pC);


  // Verts ABCD are ordered in a right-handed way.  The faces are
  // defined relative to them so that there is no possible ambiguity.
  // Note that the face definitions here DO NOT correspond to
  // iFace[ABCD](iCell).
  Vert *pVA = pC->pVVert (0);
  Vert *pVB = pC->pVVert (1);
  Vert *pVC = pC->pVVert (2);

  Face *pFA = pTC->pFFaceOpposite (pVA);
  Face *pFB = pTC->pFFaceOpposite (pVB);
  Face *pFC = pTC->pFFaceOpposite (pVC);

  // Get new cells
  Cell *pCA = pC;
  Cell *pCB = pCNewCell ();
  Cell *pCC = pCNewCell ();
  pCB->vSetRegion (pC->iRegion ());
  pCC->vSetRegion (pC->iRegion ());

  // Create new faces and set them up
  Face *pFA1 = pFNewFace ();
  Face *pFB1 = pFNewFace ();
  Face *pFC1 = pFNewFace ();
  pFA1->vSetFaceLoc (Face::eInterior);
  pFB1->vSetFaceLoc (Face::eInterior);
  pFC1->vSetFaceLoc (Face::eInterior);

  pFA1->vAssign (pCB, pCC, pVA, pVNew);
  pFB1->vAssign (pCC, pCA, pVB, pVNew);
  pFC1->vAssign (pCA, pCB, pVC, pVNew);

  // Fix old faces
  pFA->vReplaceCell (pC, pCA);
  pFB->vReplaceCell (pC, pCB);
  pFC->vReplaceCell (pC, pCC);

  // Set up new cells
  pCA->vAssign (pFA, pFB1, pFC1);
  pCB->vAssign (pFB, pFC1, pFA1);
  pCC->vAssign (pFC, pFA1, pFB1);

  assert (pCA->iFullCheck ());
  assert (pCB->iFullCheck ());
  assert (pFC->iFullCheck ());
  assert (pFA->iFullCheck ());
  assert (pFB->iFullCheck ());
  assert (pFC->iFullCheck ());
  assert (pFA1->iFullCheck ());
  assert (pFB1->iFullCheck ());
  assert (pFC1->iFullCheck ());

  // Add a vertex hint for the new vert
  pVNew->vSetHintFace (pFC1);
  pVNew->vSetType (Vert::eInterior);

  // Do some local swapping to clean up
  SUMAA_LOG_EVENT_END (INSERT_VERTEX);
  SUMAA_LOG_EVENT_BEGIN (NET_SWAPPING);
  int iRetVal = 0;
  if (qSwap) {
    iRetVal += iFaceSwap (pFA);
    iRetVal += iFaceSwap (pFB);
    iRetVal += iFaceSwap (pFC);
  }
  SUMAA_LOG_EVENT_END (NET_SWAPPING);
  if (qExitAsInsertVertex) {
    SUMAA_LOG_EVENT_BEGIN (INSERT_VERTEX);
  }
  spVUpdateLS.insert(pVNew);
  return (iRetVal);
}

//@ Add a site (or some variant if it lies outside the mesh) to the mesh
// Return value indicates success
bool Mesh2D::
qInsertPoint (const double adPoint[2], Cell * const pCGuess,
	      int *const piSwaps, const bool qSwap,
	      const bool qForce, Vert * pVNew)
{
  *piSwaps = 0;
  //@@ Find which cell the site falls within (or near, if outside the mesh)
  assert (pCGuess->iFullCheck ());
  Cell *pC = pCGuess, *pCNew = pC, *pCOld = pC;
  int iOrient = 0;
  int iNViolated, iNTies;
  Face *apFViol[2], *apFTies[2], *pF = pFInvalidFace;
  Vert VTemp;
  VTemp.vSetCoords (2, adPoint);

  if (pVNew->qValid ()) {
    assert (qForce);		// Can't twiddle the point location if it's an
    // existing vertex.

    assert (1.e-12 > dDIST2D (pVNew->adCoords (), adPoint));
  }

  do {
    int i;
    iNViolated = 0;
    iNTies = 0;
    for (i = 0; i < 3; i++) {
      pF = pC->pFFace (i);
      // If the face is pointing away from this cell, then signs need to
      // be reversed to walk in the proper direction
      int iSign = (pF->pCCellLeft () == pC) ? 1 : -1;
      iOrient = iOrient2D (pF->pVVert (0), pF->pVVert (1), &VTemp) * iSign;
      // If the vertex lies outside the current cell and behind the
      // current face, walk in that direction.  If the walk would
      // involve going outside the boundary, continue checking to try to
      // find another face to walk across; this improves behavior of
      // walks which try to insert a boundary point.
      if (iOrient == -1) {
	pCNew = pF->pCCellOpposite (pC);
	if ((pCNew->eType () == Cell::eBdryEdge) ||
	    (pF->iFaceLoc () == Face::eBdryTwoSide))
	  apFViol[iNViolated++] = pF;
	else
	  break;
      }
      if (iOrient == 0)
	apFTies[iNTies++] = pF;
      assert (iNTies >= 0 && iNTies <= 2);
      assert (iNViolated >= 0 && iNViolated <= 2);
    }
    pCOld = pC;
    pC = pCNew;
  }
  while (iOrient == -1 && pC->eType () != Cell::eBdryEdge &&
	 pF->iFaceLoc () != Face::eBdryTwoSide);

  // Restore identity of last interior cell visited, in case it's needed
  pC = pCOld;
  //@@ Check for near-coincidence with existing verts

  if (!qForce && dfLengthScale) {
    double dTestLen = dfLengthScale (adPoint) * 0.5;
    if (dDIST2D (adPoint, (pC->pVVert (0))->adCoords ()) < dTestLen ||
	dDIST2D (adPoint, (pC->pVVert (1))->adCoords ()) < dTestLen ||
	dDIST2D (adPoint, (pC->pVVert (2))->adCoords ()) < dTestLen)
      return (0);
  }

  // Define the threshold of barycentric coordinate that will induce
  // snap to face/edge
//  double dBaryThresh;
//  if (!qForce)
//    dBaryThresh = 0.16666;
//  else
//    dBaryThresh = 1.e-12;

  //@@ Decide which insertion case and insert
  switch (iNViolated) {
    //@@@ The proposed vert lies inside the mesh
  case 0:
    {
//      double adBary[3];
//      assert (pC->eType () == Cell::eTriCell);
//      ((TriCell *) pC)->vBarycentrics (adPoint, adBary);
//      assert ((adBary[0] >= -1.e-10) &&
//	      (adBary[1] >= -1.e-10) &&
//	      (adBary[2] >= -1.e-10));
//      int iLowBary = ((adBary[0] <= dBaryThresh ? 1 : 0) +
//		      (adBary[1] <= dBaryThresh ? 1 : 0) +
//		      (adBary[2] <= dBaryThresh ? 1 : 0));
//      assert (iLowBary <= 2);
      // Carl said it was fine to remove this...
//      assert (iNTies <= iLowBary);
//      switch (iLowBary) {
      switch (iNTies) {
      case 0:
	iInterior++;
	if (!(pVNew->qValid ()))
	  pVNew = pVNewVert ();
	assert (pVNew->qValid ());
	pVNew->vSetCoords (2, adPoint);
	*piSwaps = iInsertInInterior (pVNew, pC, qSwap);
	return (1);
      case 1:
	{
//	  // Face insertion
//	  if (adBary[0] <= dBaryThresh)
//	    pF = ((TriCell *) pC)->pFFaceOpposite (pC->pVVert (0));
//	  else if (adBary[1] <= dBaryThresh)
//	    pF = ((TriCell *) pC)->pFFaceOpposite (pC->pVVert (1));
//	  else if (adBary[2] <= dBaryThresh)
//	    pF = ((TriCell *) pC)->pFFaceOpposite (pC->pVVert (2));
//	  else
//	    assert2 (0, "Which barycentric was low again?");

	  // Grab the correct face for this cell
	  pF = apFTies[0];
	  if (!(pVNew->qValid ()))
	    pVNew = pVNewVert ();
	  assert (pVNew->qValid ());
	  if (qForce) {
	    pVNew->vSetCoords (2, adPoint);
	  }
	  else {
	    double adNewPoint[] =
	      {
		(pF->pVVert (0)->dX () + pF->pVVert (1)->dX ()) * 0.5,
		(pF->pVVert (0)->dY () + pF->pVVert (1)->dY ()) * 0.5
	      };
	    pVNew->vSetCoords (2, adNewPoint);
	  }
	  *piSwaps = iInsertOnFace (pVNew, pF, pC, qSwap, qForce);
	  if (*piSwaps == -1) {
	    // Tried to insert on a bdry edge when bdry's are not allowed
	    // to be changed.
	    pVNew->vMarkDeleted ();
	    *piSwaps = 0;
	    return 0;
	  }
	  return (1);
	}
      case 2:
	// Insertion where there's already a vertex.
	return (0);
      }
      assert2 (0, "Unreachable; should have returned");
      return 0;
    }
    //@@@ Outside:  Insert at edge midside
    // Only one face has the wrong orientation wrt the proposed vert
  case 1:
    {
      iOnBdryEdge++;
      pF = apFViol[0];
      double adNewPoint[] =
      {
	(pF->pVVert (0)->dX () + pF->pVVert (1)->dX ()) * 0.5,
	(pF->pVVert (0)->dY () + pF->pVVert (1)->dY ()) * 0.5
      };
      if (!(pVNew->qValid ()))
	pVNew = pVNewVert ();
      assert (pVNew->qValid ());
      pVNew->vSetCoords (2, adNewPoint);
      if (pF->iFaceLoc () == Face::eBdryFace)
	*piSwaps = iInsertOnBdryFace (pVNew, pF, pC, qSwap, qForce);
      else {
	assert (pF->iFaceLoc () == Face::eBdryTwoSide);
	*piSwaps = iInsertOnInternalBdryFace (pVNew, pF, qSwap, qForce);
      }
      return 1;
    }
    //@@@ Outside:  Logical insertion site is at existing vert; abort
  case 2:
    return (0);
  default:
    assert2 (0, "Unreachable case");
    return 0;
  }
//   // Some compilers complain without this.
//   assert2 (0, "Unreachable; should have returned");
//   return 0;
}
