C++ GEdge类代码示例

static double LC_MVertex_CURV(GEntity *ge, double U, double V)
{
  double Crv = 0;
  switch(ge->dim()){
  case 0:
    Crv = max_edge_curvature((const GVertex *)ge);
    //Crv = std::max(max_surf_curvature_vertex((const GVertex *)ge), Crv);
    // Crv = max_surf_curvature((const GVertex *)ge);
    break;
  case 1:
    {
      GEdge *ged = (GEdge *)ge;
      Crv = ged->curvature(U);
      Crv = std::max(Crv, max_surf_curvature(ged, U));
      // Crv = max_surf_curvature(ged, U);
    }
    break;
  case 2:
    {
      GFace *gf = (GFace *)ge;
      Crv = gf->curvature(SPoint2(U, V));
    }
    break;
  }
  double lc = Crv > 0 ? 2 * M_PI / Crv / CTX::instance()->mesh.minCircPoints : MAX_LC;
  return lc;
}

// compute the mesh size at a given vertex due to prescribed sizes at
// mesh vertices
static double LC_MVertex_PNTS(GEntity *ge, double U, double V)
{
  switch(ge->dim()){
  case 0:
    {
      GVertex *gv = (GVertex *)ge;
      double lc = gv->prescribedMeshSizeAtVertex();
      // FIXME we might want to remove this to make all lc treatment consistent
      if(lc >= MAX_LC) return CTX::instance()->lc / 10.;
      return lc;
    }
  case 1:
    {
      GEdge *ged = (GEdge *)ge;
      GVertex *v1 = ged->getBeginVertex();
      GVertex *v2 = ged->getEndVertex();
      if (v1 && v2){
        Range<double> range = ged->parBounds(0);
        double a = (U - range.low()) / (range.high() - range.low());
        double lc = (1 - a) * v1->prescribedMeshSizeAtVertex() +
          (a) * v2->prescribedMeshSizeAtVertex() ;

        // FIXME we might want to remove this to make all lc treatment consistent
        if(lc >= MAX_LC) return CTX::instance()->lc / 10.;
        return lc;
      }
      else
        return MAX_LC;
    }
  default:
    return MAX_LC;
  }
}

GEdgeLoop::GEdgeLoop(const std::list<GEdge*> &cwire)
{
  // Sometimes OCC puts a nasty degenerated edge in the middle of the
  // wire ...  pushing it to front fixes the problem as it concerns
  // gmsh
  std::list<GEdge*> wire;
  std::vector<GEdge*> degenerated;
  GEdge *degeneratedToInsert = 0;
  for (std::list<GEdge*>::const_iterator it = cwire.begin(); it != cwire.end(); ++it){
    GEdge *ed = *it;
    if (ed->degenerate(0))degenerated.push_back(ed);
    else wire.push_back(ed);
  }
  
  if (degenerated.size() == 1){
    wire.push_front(degenerated[0]);
  }
  else if (degenerated.size() == 2){
    degeneratedToInsert = degenerated[1];
    wire.push_front(degenerated[0]);
  }
  else if (degenerated.size() > 2){
    Msg::Error("More than two degenerated edges in one model face");
  }
  
  while (!wire.empty()){
    //    printf("wire.size = %d\n",wire.size());
    loopTheLoop(wire,loop,&degeneratedToInsert);
    //    break;
  }
}

GEdge *GNode< Point >::Search( Long n, Long i ) const {
   GEdge *p;
   for (p=adjacents; p != NULL; p=p->Next())
      if (( p->Item() == i) && (p->Node() == n))
	 break;
   return p;
}

/**
 * Destructor for node - removes edge controls on incidient edges and disconnects item from scene
 */
GNode::~GNode()
{
    graph()->invalidateRanking();
    if ( ( isEdgeControl() || isEdgeLabel())
         && isNotNullP( firstPred()) 
         && isNotNullP( firstSucc())
         && isNotNullP( firstPred()->pred())
         && isNotNullP( firstSucc()->succ()))
    {
        GRAPH_ASSERTD( areEqP( firstPred()->style(), firstSucc()->style()),
                       "Different styles on the same edge");
        GEdge *e = graph()->newEdge( firstPred()->pred(), firstSucc()->succ());
        e->setStyle( firstPred()->style());
    } else if ( isSimple())
    {
        QList< GNode *> nodes;
        GEdge* edge;
		Marker m = graph()->newMarker();
        for ( edge = firstSucc(); isNotNullP( edge); edge = edge->nextSucc())
        {
            edge->item()->adjust();
            GNode* succ = edge->succ();

            while ( succ->isEdgeControl() || succ->isEdgeLabel())
            {
                assert( isNotNullP( succ->firstSucc()));
                if ( succ->mark( m))
				{
					nodes << succ;
				}
                succ = succ->firstSucc()->succ();
            }
        }
        for ( edge = firstPred(); isNotNullP( edge); edge = edge->nextPred())
        {
            if ( edge->isSelf()) // We've already processed this one in previous loop
				continue;

			edge->item()->adjust();
            GNode* pred = edge->pred();

            while ( pred->isEdgeControl() || pred->isEdgeLabel())
            {
                assert( isNotNullP( pred->firstPred()));
                if ( pred->mark( m))
				{
					nodes << pred;
				}
                pred = pred->firstPred()->pred();
            }
        }
        
        foreach ( GNode *n, nodes)
        {
            graph()->deleteNode( n);
        }
		graph()->freeMarker( m);
    }

/**
 * Change edge's style
 */
void GGraph::setEdgeStyle( GStyle *style)
{
    if ( sel_edges.isEmpty())
        return;

    GEdge* edge = sel_edges.first();
    edge->setStyle( style);
    edge->item()->adjust();
}

GEdge * Graph::getEdge(long firstNode, long secondNode){
   GEdge *p = tnode[firstNode]->getIncidentEdges();
   while(p){
      if (p->IncidentNode() == secondNode)
	 return p;
      else
	 p=p->Next();
   }
   return NULL;      
}

bool reparamMeshVertexOnFace(MVertex *v, const GFace *gf, SPoint2 &param,
                             bool onSurface)
{

  if (gf->geomType() == GEntity::CompoundSurface ){
    GFaceCompound *gfc = (GFaceCompound*) gf;
    param = gfc->parFromVertex(v);
    return true;
  }

  if(v->onWhat()->geomType() == GEntity::DiscreteCurve ||
     v->onWhat()->geomType() == GEntity::BoundaryLayerCurve){
     param = gf->parFromPoint(SPoint3(v->x(), v->y(), v->z()), onSurface);
    return true;
  }

  if(v->onWhat()->dim() == 0){
    GVertex *gv = (GVertex*)v->onWhat();
    // hack for bug in periodic curves
    if (gv->getNativeType() == GEntity::GmshModel && gf->geomType() == GEntity::Plane)
      param = gf->parFromPoint(SPoint3(v->x(), v->y(), v->z()), onSurface);
    else
      param = gv->reparamOnFace(gf, 1);
    // shout, we could be on a seam
    std::list<GEdge*> ed = gv->edges();
    for(std::list<GEdge*>::iterator it = ed.begin(); it != ed.end(); it++)
      if((*it)->isSeam(gf)) return false;
  }
  else if(v->onWhat()->dim() == 1){
    GEdge *ge = (GEdge*)v->onWhat();
    double t;
    v->getParameter(0, t);
    param = ge->reparamOnFace(gf, t, 1);
    if(!v->getParameter(0,t)) {
      Msg::Error("Vertex %p not MEdgeVertex", v);
      return false;
      //param = gf->parFromPoint(SPoint3(v->x(), v->y(), v->z()), onSurface);
    }

    // shout, we are on a seam
    if(ge->isSeam(gf))
      return false;
  }
  else{
    double uu, vv;
    if(v->onWhat() == gf && v->getParameter(0, uu) && v->getParameter(1, vv)){
       param = SPoint2(uu, vv);
    }
    else {
      // brute force!
      param = gf->parFromPoint(SPoint3(v->x(), v->y(), v->z()), onSurface);
    }
  }
  return true;
}

double distanceToGeometry(GModel *gm, int dim, int tag, int distType,
                          double tol, int meshDiscr, int geomDiscr)
{
  double maxDist = 0.;

  if (dim == 2) {
    GEdge *ge = gm->getEdgeByTag(tag);
    if (ge->geomType() == GEntity::Line) return 0.;
    for (unsigned int i = 0; i < ge->lines.size(); i++) {
      double dist;
      switch (distType) {
        case CADDIST_TAYLOR:
          dist = taylorDistanceEdge(ge->lines[i], ge);
          break;
        case CADDIST_FRECHET:
          dist = discreteFrechetDistanceEdge(ge->lines[i], ge,
                                             tol, meshDiscr, geomDiscr);
          break;
        case CADDIST_HAUSFAST:
          dist = discreteHausdorffDistanceFastEdge(ge->lines[i], ge,
                                                   tol, meshDiscr, geomDiscr);
          break;
        case CADDIST_HAUSBRUTE:
          dist = discreteHausdorffDistanceBruteEdge(ge->lines[i], ge,
                                                    tol, meshDiscr, geomDiscr);
          break;
        default:
          Msg::Error("Wrong CAD distance type in distanceToGeometry");
          return -1.;
          break;
      }
      maxDist = std::max(dist, maxDist);
    }
  }
  else if (dim == 3) {
    if (distType == CADDIST_TAYLOR) {
      GFace *gf = gm->getFaceByTag(tag);
      if (gf->geomType() == GEntity::Plane) return 0.;
      for (unsigned int i = 0; i < gf->triangles.size(); i++)
        maxDist = std::max(taylorDistanceFace(gf->triangles[i], gf), maxDist);
      for (unsigned int i = 0; i < gf->quadrangles.size(); i++)
        maxDist = std::max(taylorDistanceFace(gf->quadrangles[i], gf), maxDist);
    }
    else {
      Msg::Error("CAD distance type %i not implemented for surfaces", distType);
      return -1.;
    }
  }
  else {
    Msg::Error("CAD distance cannot be computed for dimension %i", dim);
    return -1.;
  }

  return maxDist;
}

GEdgeSigned nextOne(GEdgeSigned *thisOne, std::list<GEdge*> &wire)
{
  if(!thisOne){
    GEdge *ge = *(wire.begin());
    wire.erase(wire.begin());
    return GEdgeSigned(1, ge);   
  }

  GVertex *gv = thisOne->getEndVertex();

  std::list<GEdge*> possibleChoices;

  std::list<GEdge*>::iterator it = wire.begin();
  std::list<GEdge*>::iterator ite = wire.end();
  while(it != ite){
    GEdge *ge = *it;
    GVertex *v1 = ge->getBeginVertex();
    GVertex *v2 = ge->getEndVertex();
    if(v1 == gv || v2 == gv) possibleChoices.push_back(ge);
    ++it;
  }
  it = possibleChoices.begin();
  ite = possibleChoices.end();
  while(it != ite){
    GEdge *ge = *it;
    if(countInList(possibleChoices, ge) == 2){
      wire.erase(std::remove_if(wire.begin(), wire.end(), 
                                std::bind2nd(std::equal_to<GEdge*>(), ge)), 
                 wire.end());
      wire.push_back(ge);
      GVertex *v1 = ge->getBeginVertex();
      GVertex *v2 = ge->getEndVertex();
      if(v1 == gv) return GEdgeSigned(1, ge);   
      if(v2 == gv) return GEdgeSigned(-1, ge);   
      Msg::Error("Something wrong in edge loop 1");
      thisOne->print();
    }
    ++it;
  }
  it = possibleChoices.begin();
  ite = possibleChoices.end();
  while(it != ite){
    GEdge *ge = *it;
    if(ge != thisOne->ge){
      wire.erase(std::remove_if(wire.begin(),wire.end(), 
                                std::bind2nd(std::equal_to<GEdge*>(), ge)), 
                 wire.end());
      GVertex *v1 = ge->getBeginVertex();
      GVertex *v2 = ge->getEndVertex();
      if(v1 == gv) return GEdgeSigned(1, ge);   
      if(v2 == gv) return GEdgeSigned(-1, ge);
      Msg::Error("Something wrong in edge loop 2");
      thisOne->print();
    }   
    ++it;
  }
  
  // should never end up here
  return GEdgeSigned(0, 0);
}

void Centerline::importFile(std::string fileName)
{
  current = GModel::current();
  std::vector<GFace*> currentFaces(current->firstFace(), current->lastFace());
  for (unsigned int i = 0; i < currentFaces.size(); i++){
    GFace *gf = currentFaces[i];
     if (gf->geomType() == GEntity::DiscreteSurface){
     	for(unsigned int j = 0; j < gf->triangles.size(); j++)
     	  triangles.push_back(gf->triangles[j]);
	if (is_cut){
	  gf->triangles.clear();
	  gf->deleteVertexArrays();
	  current->remove(gf);
	}
     }
  }

  if(triangles.empty()){
    Msg::Error("Current GModel has no triangles ...");
    return;
  }

  mod = new GModel();
  mod->load(fileName);
  mod->removeDuplicateMeshVertices(1.e-8);
  current->setAsCurrent();
  current->setVisibility(1);

  int maxN = 0.0;
  std::vector<GEdge*> modEdges(mod->firstEdge(), mod->lastEdge());
  MVertex *vin = modEdges[0]->lines[0]->getVertex(0);
  ptin = SPoint3(vin->x(), vin->y(), vin->z());
  for (unsigned int i = 0; i < modEdges.size(); i++){
    GEdge *ge = modEdges[i];
    for(unsigned int j = 0; j < ge->lines.size(); j++){
      MLine *l = ge->lines[j];
      MVertex *v0 = l->getVertex(0);
      MVertex *v1 = l->getVertex(1);
      std::map<MVertex*, int>::iterator it0 = colorp.find(v0);
      std::map<MVertex*, int>::iterator it1 = colorp.find(v1);
      if (it0 == colorp.end() || it1 == colorp.end()){
	lines.push_back(l);
	colorl.insert(std::make_pair(l, ge->tag()));
	maxN = std::max(maxN, ge->tag());
       }
      if (it0 == colorp.end()) colorp.insert(std::make_pair(v0, ge->tag()));
      if (it1 == colorp.end()) colorp.insert(std::make_pair(v1, ge->tag()));
    }
 }

  createBranches(maxN);
}

void Centerline::createSplitCompounds()
{
  //number of discrete vertices, edges, faces and regions for the mesh
  NV = current->getMaxElementaryNumber(0);
  NE = current->getMaxElementaryNumber(1);
  NF = current->getMaxElementaryNumber(2);
  NR = current->getMaxElementaryNumber(3);

  // Remesh new faces (Compound Lines and Compound Surfaces)
  Msg::Info("Centerline: creating split compounds ...");

  //Parametrize Compound Lines
  for (int i=0; i < NE; i++){
    std::vector<GEdge*>e_compound;
    GEdge *pe = current->getEdgeByTag(i+1);//current edge
    e_compound.push_back(pe);
    int num_gec = NE+i+1;
    Msg::Info("Create Compound Line (%d) = %d discrete edge",
              num_gec, pe->tag());
    GEdge *gec =  current->addCompoundEdge(e_compound,num_gec);
    if (CTX::instance()->mesh.algo2d != ALGO_2D_BAMG){
      gec->meshAttributes.method = MESH_TRANSFINITE;
      gec->meshAttributes.nbPointsTransfinite = nbPoints+1;
      gec->meshAttributes.typeTransfinite = 0;
      gec->meshAttributes.coeffTransfinite = 1.0;
    }
  }

  // Parametrize Compound surfaces
  std::list<GEdge*> U0;
  for (int i=0; i < NF; i++){
    std::vector<GFace*> f_compound;
    GFace *pf =  current->getFaceByTag(i+1);//current face
    f_compound.push_back(pf);
    int num_gfc = NF+i+1;
    Msg::Info("Create Compound Surface (%d) = %d discrete face",
              num_gfc, pf->tag());

    //1=conf_spectral 4=convex_circle, 7=conf_fe
    GFace *gfc = current->addCompoundFace(f_compound, 7, 0, num_gfc);

    gfc->meshAttributes.recombine = recombine;
    gfc->addPhysicalEntity(1);
    current->setPhysicalName("wall", 2, 1);//tag 1

  }
}

/*
 * Creates a copy of graph g.
 */
Graph3d &Graph3d::operator=( const Graph3d &g ) {
   Long	i;
   GEdge *l;

   New(g.size,g.ndep,g.nrow,g.ncol);
   for (i=0; i<size; i++){
      if (g[i] != NULL){
	 Add(i,g[i]->Item(),g[i]->seed);
	 tnode[i]->value=g[i]->value;
	 for (l=g[i]->Neighbours(); l!=NULL; l=l->Next()){
	    if (i>l->Node())
	       Link(i,l->Node(),l->Item(),l->weight);
	 }
      }
   }
   return *this;
}

static double _relocateVertex(GFace *gf, MVertex *ver,
                              const std::vector<MElement *> &lt, double tol)
{
  if(ver->onWhat()->dim() != 2) return 2.0;

  SPoint2 p1(0, 0);
  SPoint2 p2;
  if(ver->getParameter(0, p2[0])) {
    ver->getParameter(1, p2[1]);
  }
  else {
    return _relocateVertex2(gf, ver, lt, tol);
  }

  std::size_t counter = 0;
  for(std::size_t i = 0; i < lt.size(); i++) {
    for(std::size_t j = 0; j < lt[i]->getNumVertices(); j++) {
      MVertex *v = lt[i]->getVertex(j);
      SPoint2 pp;
      reparamMeshVertexOnFace(v, gf, pp);
      counter++;
      if(v->onWhat()->dim() == 1) {
        GEdge *ge = dynamic_cast<GEdge *>(v->onWhat());
        // do not take any chance
        if(ge->isSeam(gf)) return 2.0;
      }
      p1 += pp;
    }
  }
  p1 *= 1. / (double)counter;
  double worst;
  double xi = Maximize_Quality_Golden_Section(ver, gf, p1, p2, lt, tol, worst);
  // if (xi != 0) printf("xi = %g\n",xi);
  SPoint2 p = p1 * (1 - xi) + p2 * xi;
  GPoint pp = gf->point(p);
  if(!pp.succeeded()) return 2.0;
  ver->x() = pp.x();
  ver->y() = pp.y();
  ver->z() = pp.z();
  ver->setParameter(0, pp.u());
  ver->setParameter(1, pp.v());
  return worst;
}

GEdge *GNode::UnConnect( long incidentNode ) {
   GEdge *p = getIncidentEdges();
   GEdge *q;
   
   if (!p) return NULL;
   if(p->IncidentNode() == incidentNode){
      adjacents = p->Next();
      delete p;
   }
   return adjacents;
   
   while (p->Next())
      {
	 if(p->Next()->IncidentNode() == incidentNode){
	    q = p->Next();
	    p->Next(q->Next());
	    delete q;
	 }else
	    p = p->Next();
      }
   return adjacents;
}