// -> ../ut/ut_face.c





#ifdef _MSC_VER
#include "../xa/MS_Def0.h"
#endif

#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <stdarg.h>

#include <GL/gl.h>                     // GL_TRIANGLES



#include "../ut/ut_geo.h"              // Point ...
#include "../ut/ut_memTab.h"           // MemTab_..
#include "../ut/ut_itmsh.h"            // MSHIG_EDGLN_.. typedef_MemTab.. Fac3
#include "../ut/ut_face.h"             // UFA
#include "../ut/ut_txt.h"              // fnam_del
#include "../ut/func_types.h"               // UI_Func... SYM_..
#include "../ut/ut_TX.h"               // TX_Print
#include "../ut/ut_os.h"               // OS_ ..


/*

UFA_if_find_ipIn     find face covering point ipt

UFA_ipse_fac_esn     get startPoint and endPt of edge from face, edgeSeqNr

UFA_add_fac_st       add face,stat to MemTab(Fac3)
UFA_5fac_3fac_pt     make 5 faces from 3 faces and point
UFA_4fac_2fac_pt     make 4 faces from 2 faces and point
UFA_3fac_fac_pt      make 3 faces from face and point
UFA_2fac_fac_pt      make 2 faces from face and point on edge
UFA_1fac_fac_pt      modify face 

UFA_if_ck_pt_nxt     get next face using point
UFA_ifev_find_2ip    find neighbour-edge, get value from 2 points
UFA_2if_2pt          find 2 faces from 2 points opposite to common edge
UFA_ipn_chg_ip       change all points x -> y
UFA_ipn_chg_ip       change ip1 -> ip2 in a single face

UFA_esn_edgo         get esn of oriented-edge (2 ipt's)
UFA_2esn_comm_ck     check if esn1-esn2 are common edges

UFA_fnb_if_update    update face-links (keep breaklines)
UFA_fnb_set_edg      set link for edge and its reverse-edge
UFA_fnb_upd_1        update single fnb
UFA_fnb_move         update fnb for moved face

UFA_view__           display
UFA_disp_fb1         display boundary of indexed-triangle
*/



/*
//================================================================
  int UFA_fnb_set_edg (Fac3 *fnb, int ifc1, int esn1, int edgTyp) {
//================================================================
// set edge ifc1-esn1 and the reverse-edge to edgTyp
// see UFA_fnb_set_if_es UFA_fnb_set_val


  int    ifc2, esn2;


  printf("UFA_fnb_set_edg ifc=%d esn=%d val=%d\n",ifc, esn, edgTyp);

  ifc2 = 
  esn2 =
  
  if(esn1 == 1)       fnb[ifc1].i1 = edgTyp;
  else if(esn1 == 2)  fnb[ifc1].i2 = edgTyp;
  else if(esn1 == 3)  fnb[ifc1].i3 = edgTyp;



    // UFA_fnb_dump_1 (&fa[ifc], &fnb[ifc], " _setinb-fa[%d] = ",ifc);

  return 0;

}
*/



//================================================================
  int UFA_fnb_move (int ifOld, int ifNew, Fac3 *fac, Fac3 *fnb) {
//================================================================
// UFA_fnb_move         update fnb for moved face
// use BEFORE copying fac

  Fac3   *fnb1;

  fnb1 = &fnb[ifOld];

  // printf("UFA_fnb_move %d -> %d f = %d %d %d\n",ifOld,ifNew,
          // fnb1->i1,fnb1->i2,fnb1->i3);


  // change fac[fac[ifOld].i1]
  if(fnb1->i1 < 0) goto L_i2;
  UFA_ip_chg_ip (ifOld, ifNew, &fnb[fnb1->i1]);


  // change fac[fac[ifOld]].i2
  L_i2:
  if(fnb1->i2 < 0) goto L_i3;
  UFA_ip_chg_ip (ifOld, ifNew, &fnb[fnb1->i2]);



  // change fac[fac[ifOld]].i2
  L_i3:
  if(fnb1->i3 < 0) goto L_exit;
  UFA_ip_chg_ip (ifOld, ifNew, &fnb[fnb1->i3]);


  L_exit:
  return 0;

}


//================================================================
  int UFA_ip_chg_ip (int i1, int i2, Fac3 *fa) {
//================================================================
// UFA_ipn_chg_ip       change ip1 -> ip2 in a single face

  // printf("UFA_ip_chg_ip %d to %d fa = %d %d %d \n",i1,i2,fa->i1,fa->i2,fa->i3);


    if(fa->i1 == i1) fa->i1 = i2;
    if(fa->i2 == i1) fa->i2 = i2;
    if(fa->i3 == i1) fa->i3 = i2;


  return 0;

}


//================================================================
  int UFA_ipn_chg_ip (int ip1, int ip2, Fac3 *fa, int fNr) {
//================================================================
// UFA_ipn_chg_ip       change all points ip1 -> ip2

  int    i1;


  for(i1=0; i1<fNr; ++i1) {
    if(fa[i1].i1 == ip1) fa[i1].i1 = ip2;
    if(fa[i1].i2 == ip1) fa[i1].i2 = ip2;
    if(fa[i1].i3 == ip1) fa[i1].i3 = ip2;
  }


  return 0;

}

 
/*
//========================================================================
  int UFA_if_find_ipIn (int *io, int *ie, Point *pt1, Point *pa,
                        Fac3 *fa, int fNr, double *tol) {
//========================================================================
// find face covering point ipt
//
/// 
/// Output:
///   io      rc=0: faceNr; rc=1: pointIndex
///   pto     rc=0: pti on face (Z fixed); rc=1: NULL
///   retCod  0 OK; inside_face, on_face_edge:  io=faceNr, pto=point
///           1 OK; pti is on-edge; io=faceNr, ie=edgeNr, pto=point
///           2 OK; pti is on-point (pt1 == facePoint[ie]); io=ptInd; ie=ptSeqNr
///          -1 point is not inside mesh
/// 
// see also UFA_if_find_ptmsh MSH_pt_prjptmsh1
// TODO: use UFA_if_find_ipIn in MSH_pt_prjptmsh1



  int     irc, i1, ii1, ii2, ii3;
  Point2  p2i;


  p2i = UT2D_pt_pt3 (pt1);

  printf("UFA_if_find_ipIn %f %f\n",p2i.x,p2i.y);




  // check if point is in Triangle or on its boundary
  // loop tru triangles;
  for(i1=0; i1<fNr; ++i1) {
    *io = i1;
    ii1 = fa[i1].i1;
    ii2 = fa[i1].i2;
    ii3 = fa[i1].i3;
      // printf(" nxt-f%d - %d %d %d \n",i1,ii1,ii2,ii3);
    // test if point ip1 is inside Face ii1-ii2-ii3
    irc = UT2D_ck_pt_in_tria_tol ((Point2*)&pa[ii1],
                                  (Point2*)&pa[ii2],
                                  (Point2*)&pa[ii3], &p2i, tol);
      // printf(" ex pt_in_tria %d\n",irc);

    if(irc > 0) continue;    // outside ..

    //----------------------------------------------------------------
    if(irc == 0) {
      // retCod  0 OK; inside_face;
      *io = i1;
      return 0;
    }


    //----------------------------------------------------------------
    // irc: -1=on i1-i2, -2=on i2-i3, -3=on i3-i1
    if(irc >= -3) {
      // irc 1 - pti is on_edge; io=faceNr, ie=edgeNr, pto=point.
      *io = i1;
      if(irc == -1) *ie = 1;
      if(irc == -2) *ie = 2;
      if(irc == -3) *ie = 3;
      return 1;
    }


    //----------------------------------------------------------------
    // irc: -4: pi == p1; -5: pi == p2; -6: pi == p3.
    // irc 2 - pti_is_identical_with_point: io=pointIndex, pto=point.
    if(irc == -4) { *io = ii1; *ie = 1; }
    if(irc == -5) { *io = ii2; *ie = 2; }
    if(irc == -6) { *io = ii3; *ie = 3; }
    return 2;
  }


  //----------------------------------------------------------------
  // point is ouside mesh
    // printf("MSH_pt_prjptmsh1 pt is outside ..\n");
  return -1;

}
*/


//================================================================
  void UFA_ipse_fac_esn (int *ips, int *ipe, Fac3 *fac, int psn) {
//================================================================
// UFA_ipse_fac_esn     get startPoint and endPt of edge from face, edgeSeqNr
// Input:
//   psn    PointSeqNr:     1,2,3;

  // int   ii;


  if(psn == 1)       {*ips = fac->i1; *ipe = fac->i2;}
  else if(psn == 2)  {*ips = fac->i2; *ipe = fac->i3;}
  else if(psn == 3)  {*ips = fac->i3; *ipe = fac->i1;}

  // return ii;

}


//================================================================
  int UFA_fnb_upd_1 (int ifac, Fac3 *fa, Fac3 *fnb, int fNr) {
//================================================================
// UFA_fnb_upd_1        update single fnb
// write fnb[ifac].i2,i2,i3 and the fnb of its 3 neighbour-edges
// 
// see also UFA_fnb_init_1


  int    ifn, ien, ivn;

  printf("UFA_fnb_upd_1 ifac=%d fNr=%d\n",ifac,fNr);


  //----------------------------------------------------------------
  // fix i1 (esn1; 
  // get ifn=neighbour-face, ien=esn, ivn=value
  ifn = UFA_ifev_find_2ip (&ien, &ivn, fa[ifac].i1, fa[ifac].i2, fa, fnb, fNr);
  if(ifn < 0) {
    // no face found = inner-boundary (faces already deleted);
    // set this edge to inner-boundary (-3)
    // UFA_fnb_set_val (&fnb[ifac], 1, -3); // fnb[ivn].<ien> = ifac;
    UFA_fnb_set_val (&fnb[ifac], 1, -1); // fnb[ivn].<ien> = ifac;

/*
  if(fnb[ifn].st) {
    // neighbour-face ist NOT valid
    fnb[ifac].i1 = ifn; 
    UFA_fnb_set_val (&fnb[ivn], ien, ifac); // fnb[ivn].<ien> = ifac;
*/

  } else {
    // neighbour-face ist valid
    if(ivn >= 0) {
      // edge is normal
      fnb[ifac].i1 = ifn; 
      UFA_fnb_set_val (&fnb[ifn], ien, ifac); // fnb[ivn].<ien> = ifac;
    } else {
      // edge is a breakline
      fnb[ifac].i1 = ivn; 
    }
  }


  //----------------------------------------------------------------
  // fix i2 (esn2
  ifn = UFA_ifev_find_2ip (&ien, &ivn, fa[ifac].i2, fa[ifac].i3, fa, fnb, fNr);
  if(ifn < 0) {
    // no face found = inner-boundary (faces already deleted);
    // set this edge to inner-boundary (-3)
    // UFA_fnb_set_val (&fnb[ifac], 2, -3); // fnb[ivn].<ien> = ifac;
    UFA_fnb_set_val (&fnb[ifac], 2, -1); // fnb[ivn].<ien> = ifac;

/*
  if(fnb[ifn].st) {
    // neighbour-face ist NOT valid
    fnb[ifac].i2 = ifn;  
*/

  } else {
    // neighbour-face ist valid
    if(ivn >= 0) {
      // edge is normal
      fnb[ifac].i2 = ifn;
      UFA_fnb_set_val (&fnb[ifn], ien, ifac); // fnb[ivn].<ien> = ifac;
    } else {
      // edge is a breakline
      fnb[ifac].i2 = ivn;
    }
  }



  //----------------------------------------------------------------
  // fix i3 (esn3
  ifn = UFA_ifev_find_2ip (&ien, &ivn, fa[ifac].i3, fa[ifac].i1, fa, fnb, fNr);
  if(ifn < 0) {
    // no face found = inner-boundary (faces already deleted, -3);
    // or outer-boundary (-1);
    // set this edge to inner-boundary (-3)
    // UFA_fnb_set_val (&fnb[ifac], 3, -3); // fnb[ivn].<ien> = ifac;
    UFA_fnb_set_val (&fnb[ifac], 3, -1); // fnb[ivn].<ien> = ifac;

/*
  if(fnb[ifn].st) {
    // neighbour-face ist NOT valid
    fnb[ifac].i3 = ifn;  
*/

  } else {
    // neighbour-face ist valid
    if(ivn >= 0) {
      // edge is normal
      fnb[ifac].i3 = ifn;
      UFA_fnb_set_val (&fnb[ifn], ien, ifac); // fnb[ivn].<ien> = ifac;
    } else {
      // edge is a breakline
      fnb[ifac].i3 = ivn;
    }
  }


  //----------------------------------------------------------------
  fnb[ifac].st = 0;

    // TESTBLOCK
    // UFA_fnb_dump_1 (&fa[ifac], &fnb[ifac], "fnb_init_1-2 f[%d] ",ifac);
    // TESTBLOCK

  return 0;

}


//================================================================
  int UFA_ifev_find_2ip (int *enb, int *vnb, int ip1, int ip2,
                         Fac3 *fa, Fac3 *fnb, int fNr) {
//================================================================
// UFA_ifev_find_2ip    find neighbour-edge, get value from 2 points
// Input:
//   ip1,ip2  find neighbour-edge between these points
// Output:
//   enb      edge-sequence-nr of neighbour-edge
//   vnb      value of neighbour-edge (index to face with points ip1-ip2)
//   retcod   index of neighbour-face; -1=not-found
// see also UFA_fnb_get_2ip

  int    i1;


  printf("UFA_ifev_find_2ip ip1=%d ip2=%d fNr=%d\n",ip1,ip2,fNr);


  for(i1=0; i1<fNr; ++i1) {
    // test i1-i2
    if(fa[i1].i1 == ip2) {
      if(fa[i1].i2 == ip1) {*enb = 1; *vnb = fnb[i1].i1; goto L_exit;}

    // test i2-i3
    } else if(fa[i1].i2 == ip2) {
      if(fa[i1].i3 == ip1) {*enb = 2; *vnb = fnb[i1].i2; goto L_exit;}

    // test i3-i1
    } else if(fa[i1].i3 == ip2) {
      if(fa[i1].i1 == ip1) {*enb = 3; *vnb = fnb[i1].i3; goto L_exit;}
    }
  }

  i1 = -1;


  L_exit:
    printf("ex UFA_ifev_find_2ip f=%d e=%d v=%d\n",i1,*enb,*vnb);
  return i1;


}


//================================================================
  int UFA_esn_edgo (Fac3 *fc1, int ip1, int ip2) {
//================================================================
// UFA_esn_edgo         get esn of oriented-edge (2 ipt's)
// retCod:  -1       face fc1 does NOT have edge ip1-ip2
//          1|2|3    edge-sequence-nr of ip1-ip2
// see also UFA_ifac_ck_edgo


  if(fc1->i1 == ip1) {
    if(fc1->i2 == ip2) return 1;
  }

  if(fc1->i2 == ip1) {
    if(fc1->i3 == ip2) return 2;
  }

  if(fc1->i3 == ip1) {
    if(fc1->i1 == ip2) return 3;
  }

  return -1;

}



/*
//=========================================================================
  int UFA_if_ck_pt_nxt (int *ifc, int *ies, int ipt, int fNr, Fac3 *fa) {
//=========================================================================
// get next face using point
// Input:
//   ifc      index first face to check
// Output:
//   ifc      index face with next occurence of point ipt
//   ies      edgeNr; the startpoint of this edge is ipt
//   retCod   0    OK
//           -1    not found
//   
// see also UFA_nifac_ck_pt

  printf("UFA_if_ck_pt_nxt ifc=%d ipx=%d fNr=%d\n",*ifc,ipt,fNr);


  L_nxt:
    if(*ifc >= fNr) return -1;
    if(fa[*ifc].i1 == ipt) {*ies = 1; goto L_add;}
    if(fa[*ifc].i2 == ipt) {*ies = 2; goto L_add;}
    if(fa[*ifc].i3 == ipt) {*ies = 3; goto L_add;}
  
    *ifc += 1;
    goto L_nxt;


  L_add:
    printf("ex UFA_if_ck_pt_nxt %d %d\n",*ifc,*ies);

    return 0;

}


//================================================================
  int UFA_fnb_if_update (Fac3 *fa, Fac3 *fnb, int fNr) {
//================================================================
// update fnb;
// fnb-indices are invalid; negative values (breaklines) are correct;

// see UFA_fnb_init_1 UFA_fnb_BL_reset

  int      i1, ie;


  printf("UFA_fnb_if_update %d\n",fNr);
  // UFA_fnb_dump__ (fa, fnb, fNr, "in-UFA_fnb_if_update");


  // loop tru all faces
  for(i1=0; i1<fNr; ++i1) {

    // fix i1
    if(fnb[i1].i1 < 0) goto L_I2;  // negative values (breaklines) are correct
    // set neighbourface.i1 = edge i2-i1
    UFA_ifac_ck_edgo (&fnb[i1].i1, &ie, fa[i1].i2, fa[i1].i1, fa, fNr);

    // fix i2
    L_I2:
    if(fnb[i1].i2 < 0) goto L_I3;  // negative values (breaklines) are correct
    UFA_ifac_ck_edgo (&fnb[i1].i2, &ie, fa[i1].i3, fa[i1].i2, fa, fNr);

    // fix i3
    L_I3:
    if(fnb[i1].i3 < 0) continue;   // negative values (breaklines) are correct
    UFA_ifac_ck_edgo (&fnb[i1].i3, &ie, fa[i1].i1, fa[i1].i3, fa, fNr);
  }


    // UFA_fnb_dump__ (fa, fnb, fNr, "ex-UFA_fnb_if_update");
    // UFA_nfb_ck__ (fa, fnb, fNr, "ex-UFA_fnb_if_update");

  return 0;


  // if(irc) {
    // UFA_fnb_dump__ (fa, fnb, fNr, "nfb_ck");
    // return -99;
  // }

}
*/

 
//=========================================================================
  int UFA_add_fac_st (MemTab(Fac3) *fTab, int i1, int i2, int i3, int st) {
//=========================================================================
// UFA_add_fac_st       add face,stat to MemTab(Fac3)

  int       irc;
  long      l1;
  Fac3      f1;


  // if(fTab->rNr == 134)
  // printf("UFA_add_fac_st f%d = %d %d %d\n",fTab->rNr,i1,i2,i3);


  f1.i1 = i1;
  f1.i2 = i2;
  f1.i3 = i3;

  f1.st = st;

  return MemTab_sav ((MemTab*)fTab, &l1, &f1, 1);

}


//================================================================
  int UFA_disp_fb1 (int ii, Fac3 *fa, Point *pa) {
//================================================================
// display boundary of indexed-triangle 
// replacing MSH_test_disp_fb1
// TODO: UFA_disp_fb1 -> UFA_disp_opts ?         see UPAT_ipatch_disp_opts
// see GL_Disp_ipatch GL_Disp_nifac GL_Disp_cv 

  int     iCol=9;  //8=green
  Point   *p1, *p2, *p3;


  p1 = &pa[fa[ii].i1];
  p2 = &pa[fa[ii].i2];
  p3 = &pa[fa[ii].i3];

  GR_Disp_ln1 (p1, p2, iCol);
  GR_Disp_ln1 (p2, p3, iCol);
  GR_Disp_ln1 (p3, p1, iCol);

  return 0;

}


//===============================================================================
  int UFA_1fac_fac_pt (MemTab(Fac3) *fTab, int ipt, int ifc, int ie) {
//===============================================================================
// modify Face tru point ipt opposit to edge ie
// ipt lies opposit to Edge ie of face ifc
// all faces CCW.
// After return the neighbours of ie-nxt and ie-prev must be corrected.
// Input:
//   ifc      index of inputface
//   ie       esn of edge opposit to ipt 
//   ipt      point-index of new point
//
//             x
//           /  \
//          /    \
//         / ipt  \
//        /  /  \  \
//       / /      \ \
//      //   ifc    \\
//     x--------------x
//             ie


  int     ipn;



  printf("UFA_1fac_fac_pt ipt=%d ifc=%d ie=%d\n",ipt,ifc,ie);


  // get pointNr
  ipn = UFA_psn_opp_esn (ie);
    printf(" ipn=%d\n",ipn);

  // change point
  UFA_facpt_chg_psn (&fTab->data[ifc], ipt, ipn);


    // TESTBLOCK
    // printf("========= UFA_1fac_fac_pt MOD F%d \n",ifc);
    // UFA_fac_ck_sr_1 (ifc, &fTab->data[ifc], msh_pa);
    // TESTBLOCK


  return 0;

}


//===============================================================================
  int UFA_2fac_fac_pt (MemTab(Fac3) *fTab, Fac3 *fnb, int ipt, int ifc, int ie) {
//===============================================================================
// create 2 new Faces tru point ipt instead of face ifc
// ipt lies on Edge ie of face ifc
// all faces CCW.
// After return the neighbours of the new-face at e1 and e2 must be corrected.
// Input:
//   ifc      index of inputface
//   ie       esn of edge with ipt 
//   ipt      point-index of new point on edge ie on face fTab->data[ifc]
//
//              /| i3
//            / /|
//          /  / |               ie = 1                OUT: ifc-e1, fNr-e1
//       /ifc /  |
//    /   e  /   |
//   -------x-----
// i1       ipt   i2
//  i1-ipt-i3   <<< statt i1-i2-i3    <<<  ie = 1
//  ipt-i2-i3
//
//
//              /| i3
//            /  |
//          / __/ipt              ie = 2             OUT: ifc-e2, fNr-e2
//       /_ _/  e|
//    / _/   ifc |
//   -------------
// i1             i2
//  i1-i2-ipt   <<< statt i1-i2-i3    <<<  ie = 2
//  ipt-i3-i1
//
//
//              /| i3
//            /e |
//       ipt/ ifc|               ie = 3                OUT: ifc-e3, fNr-e3
//       /   \ _ |
//    /         \|
//   -------------
// i1             i2
//  ipt-i2-i3   <<< statt i1-i2-i3    <<<  ie = 3
//  ipt-i1-i2
//


  int   i1, i2, i3, fNr, oldNb, igb;
  Fac3  *f1;


  // TESTBLOCK
  // if(ipt == 21) {
  // printf("UFA_2fac_fac_pt ipt=%d ifc=%d ie=%d\n",ipt,ifc,ie);
  // UFA_fnb_dump__ (fTab->data, fnb, fTab->rNr, "2fac");
  // }
  // UFA_fnb_dump_1 (&fTab->data[ifc], &fnb[ifc], " 2fac %d ",ifc);
  // TESTBLOCK


  f1 = &fTab->data[ifc];
  i1 = f1->i1;
  i2 = f1->i2;
  i3 = f1->i3;
  igb = f1->st;
  fNr = fTab->rNr;   // this is the faceIndex of the new face
    // printf(" 2fac ifc=%d newFac=%d\n",ifc,fNr);


  if(ie == 1) {
    // make i1, ipt, i3
      // printf(" 2fac_facpt-1 f%d %d %d %d\n",ifc, ipt, i2, i3);
    // change f[ifc]
    UFA_facpt_chg_st (fTab->data, ifc, i1, ipt, i3, igb);

    // create face  ipt,i2,i3
      // printf(" 2fac_facpt-1 f%d %d %d %d\n",fTab->rNr, ipt, i2, i3);
    UFA_add_fac_st (fTab, ipt, i2, i3, igb);

    // update neighbour of new face fnb[fNr];
    fnb[fNr].i1 = fnb[ifc].i1;
    fnb[fNr].i2 = fnb[ifc].i2;
    fnb[fNr].i3 = ifc;
    fnb[fNr].st = 0;

    // update fnb[ifc]; keep i1,i3.
    oldNb = fnb[ifc].i2;  // printf(" oldNb2=%d\n",oldNb2);
    fnb[ifc].i2 = fNr;

    // update neighbour of old e2 (change ifc -> fNr)
    UFA_fnb_mod_val (&fnb[oldNb], fNr, ifc);




  } else if(ie == 2) {
    // make i1, i2, ipt
      // printf(" 2fac_facpt-2 f%d %d %d %d\n",ifc, i1, i2, ipt);
    UFA_facpt_chg_st (fTab->data, ifc, i1, i2, ipt, igb);

    // make ipt, i3, i1
      // printf(" 2fac_facpt-2 f%d %d %d %d\n",fTab->rNr, i1, ipt, i3);
    UFA_add_fac_st (fTab, ipt, i3, i1, igb);

    // update neighbour of new face fnb[fNr];
    fnb[fNr].i1 = fnb[ifc].i2;
    fnb[fNr].i2 = fnb[ifc].i3;    // 2016-02-21  
    fnb[fNr].i3 = ifc;
    fnb[fNr].st = 0;

    // update fnb[ifc]; keep i1,i2.
    oldNb = fnb[ifc].i3;  // printf(" oldNb2=%d\n",oldNb2);
    fnb[ifc].i3 = fNr;

    // update neighbour of old e3 (change ifc -> fNr)
    UFA_fnb_mod_val (&fnb[oldNb], fNr, ifc);



  } else if(ie == 3) {
    // make i1, i2, ipt
      // printf(" 2fac_facpt-3 f%d %d %d %d\n",ifc, i1, ipt, i3);
    UFA_facpt_chg_st (fTab->data, ifc, ipt, i2, i3, igb);

    // make ipt, i1, i2
      // printf(" 2fac_facpt-3 f%d %d %d %d\n",fTab->rNr, i2, i3, ipt);
    UFA_add_fac_st (fTab, ipt, i1, i2, igb);

    // update neighbour of new face fnb[fNr];
    fnb[fNr].i1 = fnb[ifc].i3;
    fnb[fNr].i2 = fnb[ifc].i1;
    fnb[fNr].i3 = ifc;
    fnb[fNr].st = 0;

    // update fnb[ifc]; keep i1,i3.
    oldNb = fnb[ifc].i1;  // printf(" oldNb2=%d\n",oldNb2);
    fnb[ifc].i1 = fNr;

    // update neighbour of old e1 (change ifc -> fNr)
    UFA_fnb_mod_val (&fnb[oldNb], fNr, ifc);
  }


    // TESTBLOCK
#ifdef DEB
    // printf("========= UFA_2fac_fac_pt MOD F%d CREATE F%d ======\n",ifc,fNr);
#endif
    // UFA_fac_ck_sr_1 (ifc, &fTab->data[ifc], msh_pa);
    // UFA_fac_ck_sr_1 (fNr, &fTab->data[fNr], msh_pa);
    // if(ifc == 80) {
    // UFA_fnb_dump_1 (&fTab->data[ifc], &fnb[ifc], " %d ",ifc);
    // UFA_fnb_dump_1 (&fTab->data[fNr], &fnb[fNr], " %d ",fNr);
    // }
    // if(ipt == 21)
    // UFA_fnb_dump__ (fTab->data, fnb, fTab->rNr, "ex 2fac");
    // if(ipt == 21) exit(0);
    // TESTBLOCK

  return 0;

}

//======================================================================
  int UFA_3fac_fac_pt (MemTab(Fac3) *fTab, Fac3 *fnb, int ipt, int ifc) {
//======================================================================
// create 3 new Faces tru point ipt instead of face ifc
// all faces CW
//
// Input-face:
//  ifc:    i1-i2-i3
//
//            i3
//          /   \
//  ne3-en3/   |  \  ne2-en2
//        /   ipt   \
//     fNr+1 /   \ fNr\
//      / /   ifc    \  \
//    i1 --------------- i2
//
// Output-faces:
//  ifc:    i1-i2-ipt
//  fNr:    i2-i3-ipt
//  fNr+1:  i3-i1-ipt


  int   irc, i1, i2, i3, fNr, ne2, ne3, en2, en3, igb;
  Fac3  *f1;


  // printf("UFA_3fac_fac_pt ifc=%d ipt=%d\n",ifc,ipt);
  // UFA_fnb_dump__ (fTab->data, fnb, fTab->rNr, "3fac");

  f1 = &fTab->data[ifc];
  fNr = fTab->rNr;
  i1 = f1->i1;
  i2 = f1->i2;
  i3 = f1->i3;
  igb = f1->st;
    // printf(" f%d = %d %d %d; new f%d and f%d\n",ifc,i1,i2,i3,
           // fTab->rNr,fTab->rNr + 1);


  irc = UFA_ife_getNf_ife (&ne2, &en2, ifc, 2, fTab->data, fnb);
  if(irc < 0) {TX_Error("UFA_3fac_fac_pt E1"); return -99;}
  // if(irc < 0) return MSH_ERR__ (2, ifc, ne2);
    // printf(" ne2=%d en2=%d\n",ne2,en2);

  irc = UFA_ife_getNf_ife (&ne3, &en3, ifc, 3, fTab->data, fnb);
  if(irc < 0) {TX_Error("UFA_3fac_fac_pt E2"); return -99;}
  // if(irc < 0) return MSH_ERR__ (2, ifc, ne3);
    // printf(" ne3=%d en3=%d\n",ne3,en3);


  UFA_facpt_chg_st (fTab->data, ifc, i1, i2, ipt, igb);

  UFA_add_fac_st (fTab, i2, i3, ipt, igb);

  UFA_add_fac_st (fTab, i3, i1, ipt, igb);

  fnb[fNr].i1     = fnb[ifc].i2;
  fnb[fNr].i2     = fNr + 1;
  fnb[fNr].i3     = ifc;
  fnb[fNr].st     = 0;

  fnb[fNr + 1].i1 = fnb[ifc].i3;
  fnb[fNr + 1].i2 = ifc;
  fnb[fNr + 1].i3 = fNr;
  fnb[fNr + 1].st = 0;

  // ifc: keep i1
  fnb[ifc].i2     = fNr;
  fnb[ifc].i3     = fNr + 1;

  // set fnb-links in faces ne2 and ne3
  if(ne2 >= 0) UFA_fnb_set_if_es (fnb, ne2, en2, fNr);
  if(ne3 >= 0) UFA_fnb_set_if_es (fnb, ne3, en3, fNr + 1);


  // set fac.st=0 (modified)
  fTab->data[ifc].st = 0;
  fTab->data[fNr].st = 0;
  fTab->data[fNr + 1].st = 0;


    // TESTBLOCK
#ifdef DEB
    // printf("========= UFA_3fac_fac_pt Pt %d MOD F%d CREATE F%d F%d ======\n",
           // ipt,ifc,fNr,fNr+1);
#endif
    // UFA_nfb_ck_f (ifc, fTab->data, fnb);
    // UFA_nfb_ck_f (fNr, fTab->data, fnb);
    // UFA_nfb_ck_f (fNr+1, fTab->data, fnb);
    // UFA_fac_ck_sr_1 (ifc, &fTab->data[ifc], msh_pa);
    // UFA_fac_ck_sr_1 (fNr, &fTab->data[fNr], msh_pa);
    // UFA_fac_ck_sr_1 (fNr+1, &fTab->data[fNr+1], msh_p2a);
    // UFA_fnb_dump__ (fTab->data, fnb, fTab->rNr, "ex 3fac");
    // UFA_fnb_dump_1 (&fTab->data[ifc], &fnb[ifc], " %d ",ifc);
    // UFA_fnb_dump_1 (&fTab->data[fNr], &fnb[fNr], " %d ",fNr);
    // UFA_fnb_dump_1 (&fTab->data[fNr+1], &fnb[fNr+1], " %d ",fNr+1);
    // TESTBLOCK


  return 0;

}


//=======================================================================
  int UFA_4fac_2fac_pt (MemTab(Fac3) *fTab, Fac3 *fnb, int ipt,
                       int f00, int e00) {
//=======================================================================
// create 4 Faces tru point ipt from face with point ipt on edge e00
// keep edge e00 of face f00; change its opposite point to ipx.
// update both neighbours (create 2 new faces)
//
//            X                    X
//          / | \                / | \
//         /  |  \              /  |  \
//    fn0 /   |   \        fn0 /   |   \
//       /    |    \       en0/ f01|f10 \
//      / f00 | f10 \        /     |     \
//      x  e00| e10  x       x ---ipt-----x
//       \   ^|     /         \ f00| f11 /
//         \  |   /             \  |   /
//           \| /   fn1           \| /   fn1 
//            X                    X     en1
//              
// Input:
//   ipt     point not yet meshed; the new center of 4 faces
//   f00     face with (common) edge e00
//   e00     the edge where ipt is on; cut also the neighbourface of this edge
// Output:
//   2 modified faces (f00,f10), 2 new faces (last 2 faces of fTab)
//   retCod  0-fNr  neighbour-face of f00 (f10)
//          -1  internal Error
//
// TODO: makes problems with edge=-2 (BL)

  int    irc, f01, f10, e10, f11;
  int    fn0, en0, fn1, en1, iex;


  // if(ipt == 27) {
  // UFA_fnb_dump__ (fTab->data, fnb, fTab->rNr, "4fac"); }
  // printf("UFA_4fac_2fac_pt ipt=%d f00=%d e00=%d\n",ipt,f00,e00);
  // UFA_fnb_dump_1 (&fTab->data[f00], &fnb[f00], "  f00=%d ",f00);


  // get f10,e10 = neighbour-face/edge of f00,e00
  // f10 = UFA_if_getNf_ife (f00, e00, fTab->data);
  irc = UFA_ife_getNf_ife (&f10, &e10, f00, e00, fTab->data, fnb);
  if(irc < 0) {TX_Error("UFA_4fac_2fac_pt E1"); return -99;}
  // if(irc < 0) return MSH_ERR__ (2, f00, f10);
    // printf(" f10=%d e10=%d\n",f10,e10);


  // new face from f00-e00
  f01 = fTab->rNr;

  // get fn0,en0 = neigbour of f00-next_edge_of_e00
  iex = UFA_esn_nxt (e00);
  irc = UFA_ife_getNf_ife (&fn0, &en0, f00, iex, fTab->data, fnb);
  if(irc < 0) {TX_Error("UFA_4fac_2fac_pt E2"); return -99;}
  // if(irc < 0) return MSH_ERR__ (2, f00, fn0);
    // printf(" fn0=%d en0=%d f01=%d\n",fn0, en0, f01);


  // mod f00, create f01
  UFA_2fac_fac_pt (fTab, fnb, ipt, f00, e00);
  // undefined neighbours at edges f00-e00 and f1-e00

  // upd. fnb[f1] prev.e00
  if(fn0 >= 0) {
    UFA_fnb_set_if_es (fnb, fn0, en0, f01);
    UFA_fnb_set_if_es (fnb, f01, 2, fn0);
  }


  if(f10 < 0) goto L_exit;


    // TESTBLOCK
    // if(ipt == 28)
    // UFA_fnb_dump__ (fTab->data, fnb, fTab->rNr, "4fac");
    // UFA_fnb_dump_1 (&fTab->data[f10], &fnb[f10], "4fac-nf %d ",f10);
    // return -99;
    // TESTBLOCK


  // new face from f10-e10
  f11 = fTab->rNr;

  // get fn1,en1 = neigbour of f10-next_edge_of_e10
  iex = UFA_esn_nxt (e10);
  irc = UFA_ife_getNf_ife (&fn1, &en1, f10, iex, fTab->data, fnb);
  if(irc < 0) {TX_Error("UFA_4fac_2fac_pt E3"); return -99;}
  // if(irc < 0) MSH_ERR__ (2, f10, fn1);
    // printf(" fn1=%d en1=%d f11=%d\n",fn1, en1, f11);


  // mod f10, create f11
  UFA_2fac_fac_pt (fTab, fnb, ipt, f10, e10);
  // undefined neighbours at edges f10-e10 and f2-e10

  // upd. fnb[f2] prev.e10
  if(fn1 >= 0) {
    UFA_fnb_set_if_es (fnb, fn1, en1, f11);
    UFA_fnb_set_if_es (fnb, f11, 2, fn1);
  }


  UFA_fnb_set_if_es (fnb, f10, e10, f01);   // << ??
  UFA_fnb_set_if_es (fnb, f11, 1, f00);
  UFA_fnb_set_if_es (fnb, f01, 1, f10);
  UFA_fnb_set_if_es (fnb, f00, e00, f11);


  L_exit:

    // TESTBLOCK
    // if(ipt == 27) {
    // UFA_fnb_dump__ (fTab->data, fnb, fTab->rNr, "ex 4fac");
    // exit(0);
    // }
    // printf(" f00=%d f01=%d e00=%d\n",f00,f01,e00);
    // printf(" f10=%d f11=%d e10=%d\n",f10,f11,e10);
    // UFA_nfb_ck_f (f00, fTab->data, fnb);
    // UFA_nfb_ck_f (f01, fTab->data, fnb);
    // UFA_nfb_ck_f (f10, fTab->data, fnb);
    // UFA_nfb_ck_f (f11, fTab->data, fnb);
    // UFA_fac_ck_sr_1 (f00, &fTab->data[f00], msh_p2a);
    // UFA_fac_ck_sr_1 (f01, &fTab->data[f01], msh_p2a);
    // UFA_fac_ck_sr_1 (f10, &fTab->data[f10], msh_p2a);
    // UFA_fac_ck_sr_1 (f11, &fTab->data[f11], msh_p2a);
    // UFA_fnb_dump__ (fTab->data, fnb, fTab->rNr, "ex 4fac");
    // return -99;
    // printf(" f1=%d f2=%d\n",f1,f2);
    // if(f00 == 4) return -99;
    // END TESTBLOCK


  return f10;

}


//=======================================================================
  int UFA_5fac_3fac_pt (MemTab(Fac3) *fTab, Fac3 *fnb, int ipt,
                       int f00, int e00) {
//=======================================================================
// keep edge e00 of face f00; change its opposite point to ipx.
// update both neighbours (create 2 new faces)
//
//                   X                                X
//                 // \\                            / | \
//          n21  / /  \ \                    n21  /   |  \
//              /  /  \  \                  en21 /    |   \
//            /   /    \  \                    / f21  | f10\
//           /f20 /    \f10\                  /   >   |  >  \
//           x   /      \   x                 x------ipt-----x
//           |  /       \   |                 |  <  /  \  <  |
//           |  /  f00   \  | n11             |f20/      \f11|  n11
//           |>/    e00   \>|                 | /    f00   \ |  en11
//           X------->------x                 X------->------x
//            
//         
// Input:
//   ipt     point not yet meshed; the new center of 5 faces
//   e00     the edge to keep of face f00; connect the other 2 edges with ipt
// Output:
//   3 modified faces, 2 new faces


  int      irc, e01, e02, ex;
  int      f10, ef10, f20, ef20, f11, f21, n11, en11, n21, en21;



  printf("UFA_5fac_3fac_pt ipt=%d f00=%d e00=%d\n",ipt,f00,e00);
  // UFA_fnb_dump__ (fTab->data, fnb, fTab->rNr, "5fac_2facpt");


  // get neighbours of f00 opposit e00 (f10,ef10 & f20,ef20)
  e01 = UFA_esn_nxt (e00);
  irc = UFA_ife_getNf_ife (&f10, &ef10, f00, e01, fTab->data, fnb);
  if(irc < 0) {TX_Error("UFA_5fac_3fac_pt E1"); return -99;}
  // if(irc < 0) MSH_ERR__ (2, f00, f10);

  e02 = UFA_esn_prv (e00);
  irc = UFA_ife_getNf_ife (&f20, &ef20, f00, e02, fTab->data, fnb);
  if(irc < 0) {TX_Error("UFA_5fac_3fac_pt E2"); return -99;}
  // if(irc < 0) MSH_ERR__ (2, f00, f20);

    printf(" f00=%d e01=%d e02=%d\n", f00, e01, e02);
    printf(" f10=%d ef10=%d\n", f10, ef10);
    printf(" f20=%d ef20=%d\n", f20, ef20);


  if((f10 < 0)||(f20 < 0)) {
    printf("**** UFA_5fac_3fac_pt E001 ipt=%d %d %d\n",ipt,f10,f20);
// TODO: use single UFA_2fac_facpt
    return -1;
  }


  // get neighbours of f10, f20
  ex = UFA_esn_nxt (ef10);
  irc = UFA_ife_getNf_ife (&n11, &en11, f10, ex, fTab->data, fnb);
  if(irc < 0) {TX_Error("UFA_5fac_3fac_pt E3"); return -99;}
  // if(irc < 0) MSH_ERR__ (2, f10, n11);

  ex = UFA_esn_nxt (ef20);
  irc = UFA_ife_getNf_ife (&n21, &en21, f20, ex, fTab->data, fnb);
  if(irc < 0) {TX_Error("UFA_5fac_3fac_pt E4"); return -99;}
  // if(irc < 0) MSH_ERR__ (2, f20, n21);

    printf(" n11=%d en11=%d\n", n11, en11);
    printf(" n21=%d en21=%d\n", n21, en21);


  f11 = fTab->rNr;   // new face from f00-e00+1
  f21 = f11 + 1;     // new face from f00-e00+2
    printf(" f11=%d f21=%d\n",f11,f21);



  //----------------------------------------------------------------
  // modify f10, create f11
  UFA_2fac_fac_pt (fTab, fnb, ipt, f10, ef10);


  //----------------------------------------------------------------
  // modify f20, create f21
  UFA_2fac_fac_pt (fTab, fnb, ipt, f20, ef20);
  // undefined neighbours at edges nfc-nec and f2-nec


  //----------------------------------------------------------------
  // modify f00
  UFA_1fac_fac_pt (fTab, ipt, f00, e00);


  if(f11 >= 0) {
    UFA_fnb_set_if_es (fnb, f11, 1, f00);
    UFA_fnb_set_if_es (fnb, f11, 2, n11);
    // in f00 update neighbour f11
    UFA_fnb_set_if_es (fnb, f00, e01, f11);
    // in n11 update neighbour f11
    if(n11 >= 0) UFA_fnb_set_if_es (fnb, n11, en11, f11);
  }


  if(f21 >= 0) {
    UFA_fnb_set_if_es (fnb, f21, 1, f10);
    UFA_fnb_set_if_es (fnb, f21, 2, n21);
    // in f00 update neighbour f11
    UFA_fnb_set_if_es (fnb, f10, ef10, f21);
    // in n11 update neighbour f11
    if(n21 >= 0) UFA_fnb_set_if_es (fnb, n21, en21, f21);
  }


  L_exit:
    // TESTBLOCK
    // UFA_nfb_ck_f (f00, fTab->data, fnb);
    // UFA_nfb_ck_f (f10, fTab->data, fnb);
    // UFA_nfb_ck_f (f11, fTab->data, fnb);
    // UFA_nfb_ck_f (f20, fTab->data, fnb);
    // UFA_nfb_ck_f (f21, fTab->data, fnb);
    // UFA_fac_ck_sr_1 (f00, &fTab->data[f00], msh_pa);
    // UFA_fac_ck_sr_1 (f10, &fTab->data[f10], msh_pa);
    // UFA_fac_ck_sr_1 (f11, &fTab->data[f11], msh_pa);
    // UFA_fac_ck_sr_1 (f20, &fTab->data[f20], msh_pa);
    // UFA_fac_ck_sr_1 (f21, &fTab->data[f21], msh_pa);
    // UFA_fnb_dump__ (fTab->data, fnb, fTab->rNr, "ex 5fac");
    // return -99;
    // printf(" f1=%d f2=%d\n",f1,f2);
    // if(f00 == 4) return -99;
    // TESTBLOCK


  return 0;

}

/*
//======================================================================
  int UFA_2esn_comm_ck (int esn1, int esn2, Fac3 *fac1, Fac3 *fac2) {
//======================================================================
// UFA_2esn_comm_ck           check if esn1-esn2 are common edges
// retCod  0=yes,  1=no
//
// esn: 


  printf("UFA_2esn_comm_ck %d %d\n",esn1,esn2);
  // printf(" fac1 = %d %d %d\n",fac1->i1,fac1->i2,fac1->i3);
  // printf(" fac2 = %d %d %d\n",fac2->i1,fac2->i2,fac2->i3);


  if(UFA_ips_fac_esn(fac1,esn1) != UFA_ipe_fac_esn(fac2,esn2)) return 1;
  if(UFA_ips_fac_esn(fac2,esn2) != UFA_ipe_fac_esn(fac1,esn1)) return 1;

  return 0;

}


//===================================================================
   int UFA_2if_2pt (int *ifc1, int *ies1, int *ifc2, int *ies2,
                    int ips, int ipe, Fac3 *fa, int fNr, Fac3 *fnb) {
//===================================================================
// find 2 faces from 2 points opposite to common edge
//   UFA_2fac_flip can create edge ips-ipe
// Input:
//   ips      startPoint
//   ipe      endPoint
// Output:
//   ifc1     face with pt ips
//   ies1     common edge ifc1-ifc2 on ifc1
//   ifc2     face with pt ipe
//   ies2     common edge ifc1-ifc2 on ifc2
//   retCod   -1   cannot find connection ips -> ipe
//             0   edge ips-ipe already exists; ipe = end of ies1
//                   ipe is opposit ies2 of ifc2
//             1   OK; es1=oppos.ips on face fc1; es2=oppos.ipe on face fc2

//
//     retCod = 1             retCod = 0
//            X                 X  
//          / | \               | \
//         /  |  \              |  \
//    fc1 /   |   \ fc2         |   \
//       /    |    \            | fc1\
//      /     |     \           | es1 \
//   ips   es1|es2   ipe      ips------ipe
//      \     |     /           | es2 /
//        \   |   /             | fc2
//          \ | /               | /
//            X                  X
//                               

//   get points of opposite edge of all faces with startPoint;
//     get opposite point of all neigbour-faces of this edges
//     if this point == endPoint end


  int   irc, i1, i2, if1, if2, ie1, ie2, ip2;


  printf("UFA_2if_2pt %d %d\n",ips,ipe);

  if1 = -1;

  // get if1 = next face with ips, get ie1 = edge opposite ips
  //   get i1, i2 = points of opposite edge
  L_nxt_0:
    ++if1;
    if(if1 >= fNr) return -1;
    ie1 = 0;
    if(fa[if1].i1 == ips) {ie1 = 2; i1=fa[if1].i2; i2=fa[if1].i3;} else
    if(fa[if1].i2 == ips) {ie1 = 3; i1=fa[if1].i3; i2=fa[if1].i1;} else
    if(fa[if1].i3 == ips) {ie1 = 1; i1=fa[if1].i1; i2=fa[if1].i2;} else
    goto L_nxt_0;
      // printf(" if1=%d ie1=%d i1=%d i2=%d\n",if1,ie1,i1,i2);


  //----------------------------------------------------------------
  // test if points of ie1 = ipe; yes: get fnb, return with retCod 0.
  if(i1 == ipe) {
    // startpoint of edge = ipe
    *ifc1 = if1;
    *ies1 = UFA_esn_prv (ie1);
    // get if2 = the neigbourface of edge ie1; ie2=common edge of if2
    irc = UFA_ife_getNf_ife (ifc2, ies2, *ifc1, *ies1, fa, fnb);
      // printf(" getNf_ife irc=%d if2=%d ie2=%d\n",irc,if2,ie2);
    irc = 0;
    goto L_exit;
  }

  //----------------------------------------------------------------
  if(i2 == ipe) {
    // endpoint of edge = ipe
    *ifc2 = if1;
    *ies2 = UFA_esn_nxt (ie1);
    // get if2 = the neigbourface of edge ie1; ie2=common edge of if2
    irc = UFA_ife_getNf_ife (ifc1, ies1, *ifc2, *ies2, fa, fnb);
      // printf(" getNf_ife irc=%d if2=%d ie2=%d\n",irc,if2,ie2);
    irc = 0;
    goto L_exit;
  }



  //================================================================
  // get if2 = the neigbourface of edge ie1; ie2=common edge of if2
  irc = UFA_ife_getNf_ife (&if2, &ie2, if1, ie1, fa, fnb);
    // printf(" getNf_ife irc=%d if2=%d ie2=%d\n",irc,if2,ie2);
  if(irc) goto L_nxt_0;

  // ip2 = point opposite to the common edge
  UFA_ipOpp_fac_esn (&ip2, &i1, &fa[if2], ie2);

  // if ip2 = ipe then end
  if(ip2 != ipe) goto L_nxt_0;

  *ifc1 = if1;
  *ies1 = ie1;

  *ifc2 = if2;
  *ies2 = ie2;

  irc = 1;

  L_exit:
      printf("ex UFA_2if_2pt irc=%d ifc1=%d ies1=%d ifc2=%d ies2=%d\n",
              irc, *ifc1, *ies1, *ifc2, *ies2);
  return irc;

}
*/

// eof
