geometrie/src/CAD4FE_PolyCurve.cpp

//---------------------------------------------------------------------------


//---------------------------------------------------------------------------
// include MAGIC Headers
#include "gestionversion.h"
#include "mg_geometrie.h"
#include "mg_gestionnaire.h"
#include "mg_arete.h"
#include "lc_point.h"
#include "tpl_fonction.h"
#include "ot_algorithme_geometrique.h"
#include "ot_decalage_parametre.h"
#include "CAD4FE_MG_ARETE_ClosestPointOn.h"


//---------------------------------------------------------------------------
// include STL headers
//---------------------------------------------------------------------------
#include <map>
#include <ostream>
#include <string>
#include <vector>
#include <sstream>
//---------------------------------------------------------------------------


#pragma hdrstop

#include "CAD4FE_PolyCurve.h"

//---------------------------------------------------------------------------

#pragma package(smart_init)

#include <math.h>

#ifdef __BORLANDC__
#pragma warn -8012      
#pragma warn -8037
#endif

using namespace CAD4FE;

PolyCurve::PolyCurve()
{

}

PolyCurve::PolyCurve(MG_ARETE * __edge)
{
        t_min = (0);
        t_max = (-1);
        InsertCurve( __edge );
}

PolyCurve::PolyCurve(MG_SOMMET * __vertex)
{
    t_min = (0);
    t_max = (0);
    lst_vertices.push_back(__vertex);
    lst_length.push_back(0);
}       

double
PolyCurve::GetLength(MG_ARETE * __edge)
{
        return __edge->get_longueur(__edge->get_tmin(), __edge->get_tmax());
}

double
PolyCurve::GetLength(unsigned __index)
{
        return lst_length[__index];
}
             
double equation_longueur(MG_ARETE & __edge,double t)
{
        PolyCurve::VerifyRefEdgeT(&__edge, t);
        double dxyz[3];
        __edge.deriver(t,dxyz);
        return sqrt(dxyz[0]*dxyz[0]+dxyz[1]*dxyz[1]+dxyz[2]*dxyz[2]);
}


MG_ARETE * PolyCurve::GetRefEdge(unsigned int __index)
{
    return lst_ref_edges[__index];
}

MG_SOMMET * PolyCurve::GetRefVertex(unsigned int __index)
{
    return lst_vertices[__index];
}

unsigned PolyCurve::GetRefEdgeCount() 
{
    return lst_ref_edges.size();
}

unsigned PolyCurve::GetRefVertexCount() 
{
    return lst_vertices.size();
}


bool PolyCurve::ContainsRefEdge(MG_ARETE * __refEdge)
{
    return ( std::find(lst_ref_edges.begin(),lst_ref_edges.end(),__refEdge) != lst_ref_edges.end());
}

bool PolyCurve::ContainsRefVertex(MG_SOMMET * __v)
{
    return ( std::find(lst_vertices.begin(),lst_vertices.end(),__v) != lst_vertices.end());
}

bool PolyCurve::Contains(MG_ELEMENT_TOPOLOGIQUE * __topo)
{
    if (__topo->get_dimension() == 0)
        return ContainsRefVertex((MG_SOMMET *)__topo);
    if (__topo->get_dimension() == 1)
        return ContainsRefEdge((MG_ARETE *)__topo);
    return false;
}       

void
PolyCurve::VerifyRefEdgeT(MG_ARETE * __edge, double & __t)
{
    double t1 = __edge->get_tmin();
    double t2 = __edge->get_tmax();

    if ( __edge->get_courbe()->est_periodique() )
    {
            if (__t > t2 )
                    __t -= __edge->get_courbe()->get_periode();
            if (__t < t1 )
                    __t += __edge->get_courbe()->get_periode();
    }
    else
    {
        double tStart=(t1<t2)?t1:t2;
        double tEnd=(t1>t2)?t1:t2;
            if (__t > tEnd )
                    __t = tEnd;
            if (__t < tStart )
                    __t = tStart;
    }
}

void
PolyCurve::inverser(double & __t, double __point[3], double precision, bool __curvilinearLength)
{
        int j;
        double s;

        double edge_length = t_max;

        // Test wether the requested point is coincident with a vertex of the polyline
        for (std::vector<MG_SOMMET *>::iterator it_vertex = lst_vertices.begin();
                it_vertex != lst_vertices.end();
                it_vertex++)
        {
                double vertex_point [3];
                MG_SOMMET * vertex = *it_vertex;
                vertex->get_point()->evaluer(vertex_point);

                for (j = 0; j < 3; j++)
                        if (fabs(vertex_point[j] - __point[j]) > edge_length*precision)
                                break;

                if ( j == 3 )
                {
                        unsigned iVertex = std::distance ( lst_vertices.begin(), it_vertex );
                        s = 0;
                        for (unsigned k = 0; k < iVertex; k++)
                               s += lst_length[k];

                        __t = s;

                        return;
                }
        }

        // If the requested point is not on a vertex
        // then test the reference edges
        unsigned iClosestEdge;
        double tClosestEdge = 0;
        double min_distance = 1E99;
        MG_ARETE * edge;
        for (std::vector <MG_ARETE *>::iterator it_edge = lst_ref_edges.begin();
                it_edge != lst_ref_edges.end();
                it_edge++)
        {
                edge = *it_edge;
                double t_edge;
                //edge->inverser ( t_edge, __point );
                MG_ARETE_ClosestPointOn func(__point, edge);
                t_edge = func.Find();
                double xEdge = (t_edge-edge->get_tmin())/(edge->get_tmax()-edge->get_tmin());
                if (xEdge>1) t_edge = edge->get_tmax();
                if (xEdge<0) t_edge = edge->get_tmin();
                double point[3];
                edge->evaluer (t_edge, point );
                double distance = OT_ALGORITHME_GEOMETRIQUE::VEC3_DISTANCE_VEC3(point, __point);
                if (distance < min_distance)
                {
                        min_distance = distance;
                        iClosestEdge = std::distance (lst_ref_edges.begin(), it_edge );
                        tClosestEdge = t_edge;
                }
        }

        Parameter_RefEdgeTToS ( tClosestEdge, lst_ref_edges[iClosestEdge], &s, __curvilinearLength);

        __t = s;

        return;
}

void
PolyCurve::VerifyS(double & __s)
{
                if ( est_periodique() )
                {
                        while (__s > get_tmax() )
                                __s -= get_periode();
                }
                else
                {
                        if (__s >  get_tmax() )
                                __s = get_tmax();
                }
}

double 
PolyCurve::RefVertex_GetS(MG_SOMMET * __refVertex)
{
    int N = RefVertex_GetIndex(__refVertex);
    double s;
    s=0;
    for (int i=0;i<N;i++)
        s += GetLength(i);

    return s;
}

void
PolyCurve::RefEdge_GetT(unsigned __iEdge,
                double __s,
                double *__t, double *__dt, bool __curvilinearLength
                )
{
        double t1,t2;
        MG_SOMMET * vertex = lst_vertices[ __iEdge ];
        MG_ARETE * edge = lst_ref_edges[ __iEdge ];

        if (vertex == edge->get_cosommet1()->get_sommet())
        {
            t1=edge->get_tmin();
            t2=edge->get_tmax();
        }
        else
        {
            t1=edge->get_tmax();
            t2=edge->get_tmin();
        }
        
        if ( __curvilinearLength == false )
        {
                *__dt = (t2  - t1) / GetLength( __iEdge );
                *__t = t1  + __s * (*__dt);
        }
        else
        {
                TPL_FONCTION1<double,MG_ARETE,double> longueur_bsp(*edge, equation_longueur);
                double s = __s;
                double increment = s/t_max*(t2-t1)/40;
                if (increment < 1E-4 * (t2-t1))
                {
                    *__dt = (t2  - t1) / GetLength( __iEdge );
                    *__t = t1  + __s * (*__dt);

                    return;
                }
                longueur_bsp.integrer_jusqua_gauss_2(t1, increment, s, __t);
                *__dt = (*__t - t1) / s;
        }
}

void
PolyCurve::RefEdge_GetS(unsigned __iEdge,
                double __t,
                double *__s, bool __curvilinearLength
                )
{
        double t1,t2;
        MG_SOMMET * vertex = lst_vertices[ __iEdge ];
        MG_ARETE * edge = lst_ref_edges[ __iEdge ];
        double period=edge->get_courbe()->get_periode();

        t1=edge->get_tmin();
        t2=edge->get_tmax();

        if (edge->get_courbe()->est_periodique())
        {
            if (__t < t1 - period*1E-6)
                    __t += period;
        }

        if ( __curvilinearLength == false)
        {
            if (vertex == edge->get_cosommet1()->get_sommet())
                *__s = lst_length[__iEdge] * (__t-t1)/(t2-t1);
            else
                *__s = lst_length[__iEdge] * (__t-t2)/(t1-t2);
        }
        else
        {
            if (vertex == edge->get_cosommet1()->get_sommet())
                *__s = edge->get_longueur(t1,__t);
            else
                *__s = edge->get_longueur(__t,t2);
        }
}

void
PolyCurve::InsertCurve(MG_ARETE *__edge)
{                
        if (lst_ref_edges.size() == 0)
        {
                lst_ref_edges.push_back ( __edge );
                lst_length.push_back ( GetLength(__edge) );
                lst_vertices.push_back ( __edge->get_cosommet1()->get_sommet());
                lst_vertices.push_back ( __edge->get_cosommet2()->get_sommet());
                t_min = 0;
                t_max = lst_length[0];
                return;
        }
        
        // Get the two vertices of this edge
        MG_SOMMET * V[2];
        int iV[2], // index of V[i] in lst_vertices[]  ie. V [ i ] == lst_vertices[ iV[ i ] ]
        iiVc=-1, // index < 2 of the common vertex in V[0,1]  ie. if (iVc != -1) V [ iiVc ] == lst_vertices[ iV[ iiVc ] ]
        iiVo=-1; // index < 2 of the other  vertex in V[0,1]  ie. if (iVo != -1) V [ iiVo ] == lst_vertices[ iV[ iiVo ] ]
        V[0] = __edge->get_cosommet1()->get_sommet();
        V[1] = __edge->get_cosommet2()->get_sommet();

        // Verify wether the vertices are already referenced in
        // the polycurve
        std::vector<MG_SOMMET*>::iterator itV[2];
        for (int i=0; i<2; i++)
        {
                itV[i] = std::find (lst_vertices.begin(), lst_vertices.end(), V[i]);
                if (itV[i] != lst_vertices.end())
                {  // V[i] is already in a curve of the polycurve
                        iV[i] = std::distance(lst_vertices.begin(), itV[i]);
                        if (iiVc == -1)
                                iiVc = i;
                        else
                        {
                            if (iV[iiVc] == 0)
                            {  // case of a periodic polycurve : two vertices are common, the
                               // common vertex should be the last vertex of this polycurve (iV[iiVc]=N-1),
                               // and the other vertex the first vertex of this polycurve (iV[iiVo]=0)
                                iiVo = iiVc;
                                iiVc = i;
                            }
                            else
                            {
                                iiVo = i;
                            }
                        }

                        if ( iV[i] != 0 && iV[i]+1 != (int)lst_vertices.size() )
                        {
                                printf (" Error in PolyCurve::InsertCurve : can't insert a MG_ARETE (extremities = %lu, %lu) not adjacent to the extremities of this PolyCurve (extremities = %lu, %lu)\n", __edge->get_cosommet1()->get_sommet()->get_id(), __edge->get_cosommet2()->get_sommet()->get_id(), get_sommet1()->get_id(), get_sommet2()->get_id());
                                return;
                        }
                }
                else               // vertex V[i] not found
                {
                        iV[i] = -1;
                        iiVo = i;
                }
        }


        if (iiVc == -1)
        {
                printf (" Error in PolyCurve::InsertCurve : the curve requested to append is not connected to the existing polycurve\n");
                return;
        }

        if ( iV [iiVc] == 0 ) // the inserted curve is adjacent to the first vertex of the polycurve
        {                                          
                std::reverse(lst_ref_edges.begin(), lst_ref_edges.end());
                std::reverse(lst_vertices.begin(), lst_vertices.end());
                std::reverse(lst_length.begin(), lst_length.end());
                InsertCurve(__edge);
        }
        else if ( iV [iiVc]+1 == lst_vertices.size() )  // the inserted curve is adjacent to the last vertex of the polycurve
        {
                lst_vertices.push_back ( V [ iiVo ] );
                lst_ref_edges.push_back ( __edge );
                lst_length.push_back ( GetLength(__edge) );
        }
        else
        {
            printf("PolyCurve::InsertCurve Error : the edge %lu is neither adjacent to end vertex nor to start vertex of polycurve %lu\n", __edge->get_id(), get_id());
            printf("Polycurve Start vertex = %lu, End vertex = %lu\n", get_sommet1()->get_id(),get_sommet2()->get_id());
            printf("Edge Start vertex = %lu, End vertex = %lu\n", V[0]->get_id(), V[1]->get_id());
            return;               
        }

        // update t_max        
        t_max = 0;
        for (std::vector<double>::const_iterator it_length = lst_length.begin();
                it_length != lst_length.end();
                it_length++)
                t_max += (*it_length);

}

void
PolyCurve::Parameter_SToRefEdgeT (double __s, unsigned * __iEdge, double *__t, double *__dt, bool __curvilinearLength)
{
        double s;
        unsigned i;

        VerifyS(__s);

        i=0;
        s=0;
        while (s + GetLength(i) < __s && i+1 < lst_ref_edges.size() )
        {
                s += GetLength(i);
                i++;
        }

        RefEdge_GetT(i, __s - s, __t, __dt, __curvilinearLength);
        *__iEdge = i;
}

void
PolyCurve::Parameter_SToRefEdgeT (double __s, MG_ARETE ** __edge, double *__t, double *__dt, bool __curvilinearLength)
{
    unsigned iEdge;
    Parameter_SToRefEdgeT(__s, &iEdge, __t, __dt, __curvilinearLength);
    *__edge = GetRefEdge(iEdge);
}

unsigned 
PolyCurve::RefEdge_GetIndex(MG_ARETE * __edge)
{
    unsigned i = std::distance ( lst_ref_edges.begin(), std::find ( lst_ref_edges.begin(), lst_ref_edges.end(), __edge) );
    return i;
}

unsigned 
PolyCurve::RefVertex_GetIndex(MG_SOMMET * __vertex)
{
    unsigned i = std::distance ( lst_vertices.begin(), std::find ( lst_vertices.begin(), lst_vertices.end(), __vertex) );
    return i;
}

void
PolyCurve::Parameter_RefEdgeTToS (double __t, MG_ARETE * __edge, double *__s, bool __curvilinearLength)
{                
        double s, ds;
        unsigned i;

        i = RefEdge_GetIndex(__edge);

        s=0;
        for ( unsigned k = 0; k < i; k++)
        {
                s += GetLength(k);
        }
        RefEdge_GetS(i, __t, &ds, __curvilinearLength);
        s += ds;

        *__s = s;
}
                                         
void
PolyCurve::evaluer (double __s, double __X[3])
{
        bool curvilinearLength = false;
        evaluer(__s, __X, curvilinearLength);
}

void
PolyCurve::evaluer (double __s, double __X[3], bool __curvilinearLength) // curvilinear abcissa
{
        if (IsPoint())
            {
                lst_vertices[0]->get_point()->evaluer(__X);
                return;
            }

        VerifyS (__s);
        unsigned iEdge;
        MG_ARETE * edge;
        double t, dt;
        Parameter_SToRefEdgeT ( __s, & iEdge , &t, &dt, __curvilinearLength );
        edge = lst_ref_edges [ iEdge ]; 
        PolyCurve::VerifyRefEdgeT(edge, t); 
        edge->evaluer ( t, __X );
}
                                                                                           
void
PolyCurve::deriver (double __s, double __X[3]) // curvilinear abcissa
{
        bool curvilinearLength = false;
        deriver(__s, __X, curvilinearLength);
}

void
PolyCurve::deriver (double __s, double __X[3], bool __curvilinearLength) // curvilinear abcissa
{                       
        if (IsPoint())
            {
                for (int i=0;i<3;i++) __X[i] = 1;
                return;
            }
             
        VerifyS (__s);
        unsigned iEdge;
        MG_ARETE * edge;
        double t, dt;
        Parameter_SToRefEdgeT ( __s, & iEdge , &t, &dt, __curvilinearLength );
        edge = lst_ref_edges [ iEdge ];   
        PolyCurve::VerifyRefEdgeT(edge, t);
        edge->deriver ( t, __X );
        
        // dC/ds = dC/dt * dt/ds
        for (int i = 0; i<3; i++)
                __X[i] *= dt;
}

                                                                                                                                          
void
PolyCurve::deriver_seconde (double __s, double __ddxyz[3], double *__dxyz, double *__xyz) // curvilinear abcissa
{
        bool curvilinearLength = false;
        deriver_seconde(__s, __ddxyz, __dxyz, __xyz, curvilinearLength);

}
void
PolyCurve::deriver_seconde (double __s, double __ddxyz[3], double *__dxyz, double *__xyz, bool __curvilinearLength) // curvilinear abcissa
{                   
        unsigned iEdge;
        MG_ARETE * edge;
        double t, dt;
      
        if (IsPoint())
            {
                if (__ddxyz) for (int i=0;i<3;i++) __ddxyz[i]=1;
                if (__dxyz) for (int i=0;i<3;i++) __dxyz[i]=1;
                evaluer(0,__xyz);
                return;
            }

        VerifyS (__s);

        Parameter_SToRefEdgeT ( __s, & iEdge , &t, &dt, __curvilinearLength );
        edge = lst_ref_edges [ iEdge ];

        for (int i=0; i<3; i++)
                __ddxyz[i] = 0;

        if (__dxyz)
                edge->deriver(t, __dxyz);
        // dC/ds = dC/dt * dt/ds
        for (int i = 0; i<3; i++)
                __dxyz[i] *= dt;

        if (__xyz)
                edge->evaluer(t, __xyz);
}


void
PolyCurve::Merge( PolyCurve & __polycurve)
{
        if (lst_vertices.size() == 0)
        {
                for (std::vector < MG_ARETE * >::const_iterator it = __polycurve.lst_ref_edges.begin();
                        it != __polycurve.lst_ref_edges.end();
                        it++ )
                                InsertCurve ( *it );
        }
        else if (get_sommet2() == __polycurve.get_sommet2() || get_sommet1() == __polycurve.get_sommet2() )
        {;
                for (std::vector < MG_ARETE * >::reverse_iterator it = __polycurve.lst_ref_edges.rbegin();
                        it != __polycurve.lst_ref_edges.rend();
                        it++)
                                InsertCurve ( *it );
        }
        else if ( get_sommet1() == __polycurve.get_sommet1() || get_sommet2() == __polycurve.get_sommet1() )
        {
                for (std::vector < MG_ARETE * >::const_iterator it = __polycurve.lst_ref_edges.begin();
                        it != __polycurve.lst_ref_edges.end();
                        it++ )
                                InsertCurve ( *it );
        }
        else
        {
            printf("PolyCurve::Merge Error : the polycurve %lu is neither adjacent to end vertex nor to start vertex of polycurve %lu\n", __polycurve.get_id(), get_id());
            printf("Polycurve Start vertex = %lu, End vertex = %lu\n", get_sommet1()->get_id(),get_sommet2()->get_id());
            printf("Other __polycurve Start vertex = %lu, End vertex = %lu\n", __polycurve.get_sommet1()->get_id(), __polycurve.get_sommet2()->get_id());
            return;  

        }
}

void
PolyCurve::SetPeriodicPoleRefVertex(MG_SOMMET * __v)
{
    if ( est_periodique() == false)
        return;
        
    std::vector<MG_ARETE*> e;
    std::vector<MG_SOMMET*> v;
    std::vector<double> l;
    unsigned k=RefVertex_GetIndex(__v);
    unsigned N=lst_vertices.size();
    for (unsigned i=0; i<N; i++)
    {
        unsigned m = (i+k)%N;  
        unsigned n = (i+k+1)%N;           
        MG_SOMMET * tempv = lst_vertices[ m ];
        MG_SOMMET * lastV = (v.size()) ? v[v.size()-1] : 0;
        if ( lastV && lastV == tempv )
            tempv = lst_vertices[n];
        v.push_back(tempv);
    }
    for (unsigned i=0; i+1<N; i++)
        e.push_back(lst_ref_edges[ (i+k)%(N-1) ]);
    for (unsigned i=0; i+1<N; i++)
        l.push_back(lst_length[ (i+k)%(N-1) ]);
    lst_vertices = v;
    lst_ref_edges = e;
    lst_length = l;
}

double
PolyCurve::get_longueur (double __s_min, double __s_max, double precision) 
{
        double result;
        if ( __s_min >= 0 && __s_max >= 0 )
                return fabs (__s_min - __s_max);
        else
        {
                return t_max;
        }
}

double PolyCurve::get_tmin()
{
        return t_min;
}

double PolyCurve::get_tmax()
{
        return t_max;
}

MG_SOMMET *
PolyCurve::get_sommet1() 
{
        return lst_vertices[0];
}


MG_SOMMET *
PolyCurve::get_sommet2() 
{
        return lst_vertices[lst_vertices.size()-1];
}

double
PolyCurve::get_sommet1_s()
{
        return get_tmin();
}


double
PolyCurve::get_sommet2_s() 
{
        return get_tmax();
}

void
PolyCurve::inverser(double & __t, double __point[3], double precision)
{
        bool curvilinearLength = false;

        inverser(__t, __point, precision, curvilinearLength);
}

int PolyCurve::est_periodique(void)
{
        if ( get_sommet1 () == get_sommet2 () )
                return 1;
        else
                return 0;
}

double PolyCurve::get_periode(void)
{
    if ( est_periodique() )
        return get_longueur();
    else
        return 0;
}


int PolyCurve::get_type_geometrique(TPL_LISTE_ENTITE<double> &param)
{
        return 4343; // TYPE = Poly curve
}

void PolyCurve::enregistrer(std::ostream& o,double version)
{
    unsigned int i;

    if (IsPoint() == false)
        {
    o << "%" << get_id()
      << "=CAD4FE_POLYCURVE(" << GetRefEdgeCount() << ",(";

    for (i=0; i < GetRefEdgeCount(); i++)
    {
        if ( i ) o <<",";
      o << "$" << GetRefEdge(i)->get_id();
    }                    

    o << "));" << std::endl;
        }
    else // The polycurve is a point !
        {
          o << "%"   << get_id() << "=CAD4FE_POLYCURVE(0,($"<<GetRefVertex(0)->get_id()<<"));"<< std::endl;
        }
}

 std::vector<MG_ARETE*> &
PolyCurve::GetRefEdges() 
{
        return lst_ref_edges;
}

void
CAD4FE::SplitRefEdge(MG_ARETE * __refEdge, MG_SOMMET * __refVertex1, MG_SOMMET * __refVertex2, double __xyz[3], MG_VOLUME * __refBody, MG_GEOMETRIE * __geom, MG_ARETE * edges[2], MG_SOMMET ** __splitVertex)
{
    int i,j;
    double tSplit; 
    MG_ARETE * splitEdge = __refEdge;
    splitEdge->inverser(tSplit, __xyz);

    // Create the split vertex at position xyz
    MG_POINT * point = new LC_POINT(__xyz);
    __geom->ajouter_mg_point(point);
    std::ostringstream idoriginal;
    idoriginal << "Vertex of "<< splitEdge->get_idoriginal() <<" split at t="<<tSplit;
    MG_SOMMET * splitVertex = new MG_SOMMET(idoriginal.str(), point);
    *__splitVertex = splitVertex;
    __geom->ajouter_mg_sommet(splitVertex);

    // split the edge in two edges
    MG_COURBE * curves[2];
    MG_SOMMET * vertices[2][2];
    MG_COSOMMET * cov[2][2];
                                      
    vertices[0][0] = splitEdge->get_cosommet1()->get_sommet();
    vertices[0][1] = splitVertex;
    vertices[1][0] = splitVertex;
    vertices[1][1] = splitEdge->get_cosommet2()->get_sommet();
                                                   
    curves[0] = curves[1] = splitEdge->get_courbe();

    double verticesT[2][2];
    verticesT[0][0] = splitEdge->get_tmin(); 
    verticesT[1][1] = splitEdge->get_tmax();
    verticesT[0][1] = verticesT[1][0] = tSplit;

    //  create two edges
    for (i=0; i<2; i++)
    {
        // orientation
         int orientation = splitEdge->get_orientation();
        // Create an edge using the start vertex and the split vertex
        std::ostringstream idoriginal;
        idoriginal << __refEdge->get_idoriginal();
        idoriginal << " split with t0 " << verticesT[i][0] << " t1 " << verticesT[i][1];
        edges[i] = new MG_ARETE (idoriginal.str(), 0, curves[i], orientation);
        __geom->ajouter_mg_arete (edges[i]);
        // create covertices
        for (j=0;j<2;j++)
        {
            cov[i][j] = __geom->ajouter_mg_cosommet (edges[i], vertices[i][j]);
        }

        edges[i]->changer_cosommet1(cov[i][0]);
        edges[i]->changer_cosommet2(cov[i][1]);
    }

    MG_BOUCLE * split_edge_loops[10];
    int nb_split_edge_loops = 0;
    MG_COARETE * split_edge_coedges[10], * coedges[20];
    int nb_split_edge_coedges = 0, nb_coedges = 0;
    // create four coedges (in manifold body case)
    // replace each coedges belong to initial edge in both loops of the ref faces
    //
        // for each shell of the ref body
        unsigned nb_shells = __refBody->get_nb_mg_coquille ();
        for (unsigned it_shells = 0; it_shells < nb_shells; it_shells++)
        {
                MG_COQUILLE * shell = __refBody->get_mg_coquille(it_shells);
        // for each coface->face F of the shell
                unsigned nb_coface = shell->get_nb_mg_coface();
                for (unsigned it_coface = 0; it_coface < nb_coface; it_coface++)
                {
                        MG_FACE * face = shell->get_mg_coface(it_coface)->get_face();

        // for each loop L of this face
                        unsigned nb_loop = face->get_nb_mg_boucle();
                        for (unsigned it_loop = 0; it_loop < nb_loop; it_loop++)
                        {
                                MG_BOUCLE * loop = face->get_mg_boucle(it_loop);
        // for each coedge coe of L
                                unsigned nb_edge = loop->get_nb_mg_coarete();
                                for (unsigned it_edge = 0; it_edge < nb_edge; it_edge++)
                                {
                                        MG_COARETE * coedge = loop->get_mg_coarete(it_edge);
                                        MG_ARETE * edge = coedge->get_arete();
                                          
        // if coe contains split_edge
                                        if (edge == splitEdge)
                                        {
                                            split_edge_loops[ nb_split_edge_loops++ ] = loop;        
                                            split_edge_coedges[ nb_split_edge_coedges++ ] = coedge;
        // create coe[0] ( edges[0], F, coe->sense )
        // create coe[1] ( edges[1], F, coe->sense )
        // delete coe in L
        // add coe[0] in L
        // add coe[1] in L                                            
                                            for (i=0; i<2; i++)
                                            {
                                                coedges[nb_coedges++] = __geom->ajouter_mg_coarete (edges[i], loop, coedge->get_orientation());
                                                loop->ajouter_mg_coarete(coedges[nb_coedges-1]);
                                            }

                                            loop->supprimer_mg_coarete(coedge);
                                            __geom->supprimer_mg_coarete(coedge);
                                        }
        // fi
                                }
                        }
                }
        }                          

}

void
CAD4FE::SplitPolyCurve(PolyCurve * __polyCurve, double __xyz[3], MG_VOLUME * __refBody, MG_GEOMETRIE * __geom, PolyCurve * __result[2], MG_ARETE **__origRefEdge,  MG_SOMMET ** __splitRefVertex, MG_ARETE * __splitRefEdges[2])
{
    int i,j;
    *__splitRefVertex = NULL;
    __splitRefEdges[0]=__splitRefEdges[1]=NULL;

    // get the parameter of __splitVertex
    double sSplitVertex;
    __polyCurve->inverser(sSplitVertex,__xyz,1E-6);

    unsigned index; // index of the reference edge split
    double tEdge, dtEdge; // real parameter of the reference edge split
    __polyCurve->Parameter_SToRefEdgeT(sSplitVertex,&index,&tEdge,&dtEdge, false);
    *__origRefEdge = __polyCurve->GetRefEdge(index);
    MG_SOMMET * v[2];
    v[0] = __polyCurve->GetRefVertex(index);
    v[1] = __polyCurve->GetRefVertex(index+1);

    // Test wether __xyz is located on one extremity of the ref edge
    for (i=0; i<2; i++)
    {
        double v_xyz[3];

        v[i]->get_point()->evaluer(v_xyz);

        if (OT_ALGORITHME_GEOMETRIQUE::VEC3_DISTANCE_VEC3(__xyz, v_xyz) / __polyCurve->GetLength(index) < 1E-6)
        {
            if (index+i == 0 || index+i+1 == __polyCurve->GetRefVertexCount())
            {
                printf("Can't split a polycurve on vertex ID %lu because its extremities are vertices : %lu and %lu\n", v[i]->get_id(), __polyCurve->GetRefVertex(0)->get_id(), __polyCurve->GetRefVertex(__polyCurve->GetRefVertexCount()-1)->get_id());
                return;
            }

            *__origRefEdge = NULL;
            *__splitRefVertex = v[i];
            __splitRefEdges[0] = __polyCurve->GetRefEdge(index-1+i);
            __splitRefEdges[1] = __polyCurve->GetRefEdge(index+i);
                      
            // create the two polyCurves
            for (int k=0; k<2; k++)
            {
                // create polycurve
                __result[k] = new PolyCurve ( __splitRefEdges[k] );
                __geom->ajouter_mg_courbe(__result[k]);
                for (j = (k==0)?index-2+i:index+i+1; j >= 0 && j < __polyCurve->GetRefEdgeCount(); j += -1+2*k)
                    __result[k]->InsertCurve(__polyCurve->GetRefEdge(j));
            }
        
            break;
        }
    }
    if (*__splitRefVertex == NULL)
    {
        CAD4FE::SplitRefEdge(*__origRefEdge, v[0], v[1], __xyz, __refBody, __geom, __splitRefEdges, __splitRefVertex);

        if ( __splitRefEdges[0]->get_cosommet1()->get_sommet() != v[0] && __splitRefEdges[0]->get_cosommet2()->get_sommet() != v[0] )
            std::swap(__splitRefEdges[0], __splitRefEdges[1] );  

        // create the two polyCurves
        for (i=0; i<2; i++)
        {
            // create polycurve
            __result[i] = new PolyCurve ( __splitRefEdges[i] );
            __geom->ajouter_mg_courbe(__result[i]);
            for (j = index-1+2*i; j >= 0 && j < __polyCurve->GetRefEdgeCount(); j += -1+2*i)
                __result[i]->InsertCurve(__polyCurve->GetRefEdge(j));
        }
    }

    if (__result[0]->get_sommet1() !=__result[1]->get_sommet1() &&__result[0]->get_sommet1() !=__result[1]->get_sommet2() &&
    __result[0]->get_sommet2() !=__result[1]->get_sommet1() &&__result[0]->get_sommet2() !=__result[1]->get_sommet2())
        printf("Error");

}

bool PolyCurve::IsPoint() const
{
    return lst_vertices.size() == 1 && lst_ref_edges.size() == 0;
}

void
PolyCurve::get_param_NURBS(int& indx_premier_ptctr,TPL_LISTE_ENTITE<double> &param)
{
return;
}


Revision:	763
Committed:	Wed Dec 2 19:55:53 2015 UTC (9 years, 5 months ago) by francois
File size:	30529 byte(s)
Log Message:	Le fichier MAGiC est maintenant versionné. LA version actuelle est 2.0. L'ancienne version est 1.0. Tout est transparent pour l'utilisateur. Les vieilles versions sont lisibles mais les nouveaux enregistrements sont dans la version la plus récente. Changement des conditions aux limites : ajout d'un parametre pour dire si la condition numerique est une valeur ou une formule ou un lien vers une autre entité magic. Les parametres pour saisir sont maintenant -ccf -ccfi -ccff -ccft -ccfit -ccfft