// Gnurbs - A curve and surface library
// Copyright (C) 2008-2026 Eric Bechet
//
// See the LICENSE file for contributions and license information.
// Please report all bugs and problems to <bechet@cadxfem.org>.
//

#include <iostream>
#include <vector>

#include "nvtkdisplay.h"
#include "nbeziercurve.h"
#include "nbsplinesurface.h"

// helper function to initialize boundary curves
void init_crvs(nbeziercurve &crvU1,nbeziercurve &hodoU1,nbeziercurve &crvU2,nbeziercurve &hodoU2,
                  nbeziercurve &crvV1,nbeziercurve &hodoV1,nbeziercurve &crvV2,nbeziercurve &hodoV2);


// function that computes a bicubic coons patch.
void coons3(nbeziercurve &crvU1,nbeziercurve &hodoU1,nbeziercurve &crvU2,nbeziercurve &hodoU2,
            nbeziercurve &crvV1,nbeziercurve &hodoV1,nbeziercurve &crvV2,nbeziercurve &hodoV2,
            nbsplinesurface &surf)
{
  // to code
}

int main(void)
{
  data_container data;
  nvtkdisplay display;

  const int nCP=7;
  nbeziercurve crv_u1(nCP),crv_u2(nCP); // curves and their hodograph
  nbeziercurve crv_v1(nCP),crv_v2(nCP);
  nbeziercurve hodo_u1(nCP),hodo_u2(nCP);
  nbeziercurve hodo_v1(nCP),hodo_v2(nCP);

  nbsplinesurface surface_coons;

  init_crvs(crv_u1,hodo_u1,crv_u2,hodo_u2,crv_v1,hodo_v1,crv_v2,hodo_v2);

  coons3(crv_u1,hodo_u1,crv_u2,hodo_u2,crv_v1,hodo_v1,crv_v2,hodo_v2,surface_coons); // call to the function computing the bicubic patch


  // affichage des crvs et du patch, ainsi que les derivees sur le bord.
  crv_u1.Display(data);
  crv_u2.Display(data);
  crv_v1.Display(data);
  crv_v2.Display(data);
  surface_coons.Display(data);

  for (double u=crv_u1.min_u();u<=crv_u1.max_u()+0.005;u+=0.01)
  {
    npoint p;
    npoint dp;
    npoint3 p1;
    line l;

    crv_u1.P(u,p);
    hodo_u1.P(u,dp);
    p1=p;p1+=npoint3(dp)*-0.05;
    l.pts[0]=p;l.pts[1]=p1;
    data.add_line(l);

    crv_u2.P(u,p);
    hodo_u2.P(u,dp);
    p1=p;p1+=npoint3(dp)*+0.05;
    l.pts[0]=p;l.pts[1]=p1;
    data.add_line(l);
  }

  for (double v=crv_v1.min_u();v<=crv_v1.max_u()+0.005;v+=0.01)
  {
    npoint p;
    npoint dp;
    npoint3 p1;
    line l;

    crv_v1.P(v,p);
    hodo_v1.P(v,dp);
    p1=p;p1+=npoint3(dp)*-0.05;
    l.pts[0]=p;l.pts[1]=p1;
    data.add_line(l);

    crv_v2.P(v,p);
    hodo_v2.P(v,dp);
    p1=p;p1+=npoint3(dp)*+0.05;
    l.pts[0]=p;l.pts[1]=p1;
    data.add_line(l);
  }

  display.init_data(data);
  display.display();
  return 0;
}


void init_crvs(nbeziercurve &crvU1,nbeziercurve &hodoU1,nbeziercurve &crvU2,nbeziercurve &hodoU2,
                  nbeziercurve &crvV1,nbeziercurve &hodoV1,nbeziercurve &crvV2,nbeziercurve &hodoV2)
{
  int degree=crvU1.degree();
  // 4 positional curves
  for (int i=0;i<degree+1;++i)
  {
    crvU1.CP(i)=npoint(i,0,sin(i*1.0/(degree)*n_pi));
    crvU2.CP(i)=npoint(i,degree,-sin(i*4.0/(degree)*n_pi));
    crvV1.CP(i)=npoint(0,i,sin(i*2.0/(degree)*n_pi));
    crvV2.CP(i)=npoint(degree,i,-sin(i*3.0/(degree)*n_pi));
  }
  // les 8 corner derivatives
  npoint t;
  crvU1.dPdu(0,t);t.w()=1.0; npoint3 dcrvU1du0=t;
  crvU1.dPdu(1,t);t.w()=1.0; npoint3 dcrvU1du1=t;
  crvU2.dPdu(0,t);t.w()=1.0; npoint3 dcrvU2du0=t;
  crvU2.dPdu(1,t);t.w()=1.0; npoint3 dcrvU2du1=t;
  crvV1.dPdu(0,t);t.w()=1.0; npoint3 dcrvV1dv0=t;
  crvV1.dPdu(1,t);t.w()=1.0; npoint3 dcrvV1dv1=t;
  crvV2.dPdu(0,t);t.w()=1.0; npoint3 dcrvV2dv0=t;
  crvV2.dPdu(1,t);t.w()=1.0; npoint3 dcrvV2dv1=t;
  // 4 boundary derivatives (hodograph)
  for (int i=0;i<degree+1;++i)
  {
    hodoU1.CP(i)=dcrvV2dv0*i*1.0/degree+dcrvV1dv0*(degree-i)*1.0/degree+10*npoint3(0,0,sin(i*4.0/(degree)*n_pi));
    hodoU2.CP(i)=dcrvV2dv1*i*1.0/degree+dcrvV1dv1*(degree-i)*1.0/degree+10*npoint3(0,0,-sin(i*4.0/(degree)*n_pi));
    hodoV1.CP(i)=dcrvU2du0*i*1.0/degree+dcrvU1du0*(degree-i)*1.0/degree+10*npoint3(0,0,-sin(i*4.0/(degree)*n_pi));
    hodoV2.CP(i)=dcrvU2du1*i*1.0/degree+dcrvU1du1*(degree-i)*1.0/degree+10*npoint3(0,0,sin(i*4.0/(degree)*n_pi));
  }
}