regular_triangulation_3.h

很多二维 三维几何计算算法 C++ 类库

    } else {
      return kdel_.simulator()->template event<Non_vertex_event>(redundant_points_.find(k)->second).cell_handle();
    }
  }

  void set_listener(Listener *l) {
    listener_= l;
  }
  Listener* listener() const
  {
    return listener_;
  }
  void audit_structure() const
  {
    if (!has_certificates()) {
      for (typename RPMap::const_iterator it= redundant_points_.begin(); 
	   it != redundant_points_.end(); ++it) {
	CGAL_assertion(it->second==Event_key());
      }

      for (typename RCMap::const_iterator it= redundant_cells_.begin(); 
	   it != redundant_cells_.end(); ++it) {
	Cell_handle h= it->first;
	Point_key k=it->second;
	Cell_handle lh= triangulation().locate(k);
	if (lh != h) {
	  typename Instantaneous_kernel::Current_coordinates cco= triangulation().geom_traits().current_coordinates_object();
	  bool found=false;
	  for (unsigned int i=0; i<4; ++i) {
	    if (lh->vertex(i)->point() != Point_key() && cco(lh->vertex(i)->point()).point() 
		== cco(k).point()) {
	      found=true;
	    }
	  }
	  CGAL_assertion(found);
	}

	
      }
    } else {
      for (typename Triangulation::Finite_vertices_iterator vit= triangulation().finite_vertices_begin();
	   vit != triangulation().finite_vertices_end(); ++vit) {
	if (triangulation().degree(vit) == 4) {
	  CGAL_assertion_code(Point_key k= vit->point());
	  // it could be infinite
	  // !! for VC
	  CGAL_assertion(vit->info() != Event_key() || !k.is_valid());
	}
	else {
	  CGAL_assertion_code(Point_key k= vit->point());
	  CGAL_assertion(vit->info() == Event_key());
	}
	CGAL_assertion(redundant_points_.find(vit->point())== redundant_points_.end());
      }
      CGAL_assertion(triangulation().infinite_vertex()->info() == Event_key());
      for (typename RPMap::const_iterator it= redundant_points_.begin(); it != redundant_points_.end(); ++it) {
	CGAL_assertion(kdel_.vertex_handle(it->first) == Vertex_handle());
	Cell_handle ch= get_cell_handle(it->first);
	typename RCMap::const_iterator beg= redundant_cells_.lower_bound(ch);
	typename RCMap::const_iterator end= redundant_cells_.upper_bound(ch);
	bool found=false;
	for (; beg != end; ++beg) {
	  if (beg->second == it->first) {
	    found=true;
	    break;
	  }
	}
	CGAL_assertion(found);
      } 

      for (typename RCMap::const_iterator it= redundant_cells_.begin(); 
	   it != redundant_cells_.end(); ++it) {
	Point_key pk= it->second;
	Cell_handle ch= it->first;
	CGAL_assertion(redundant_points_.find(pk) != redundant_points_.end());
	Event_key k= redundant_points_.find(pk)->second;

	Cell_handle ech= get_cell_handle(pk);
	CGAL_assertion(ch== ech);
      }
    }
  }
protected:
  //! also much check for vertex_events
  void flip(const typename Triangulation::Edge &edge) {
    typename Triangulation::Facet f= kdel_.flip(edge);

    on_geometry_changed();
  }

  void flip(const typename KDel::Facet &flip_facet) {
    typename Triangulation::Edge edge=  kdel_.flip(flip_facet);

    on_geometry_changed();
  }

  CGAL::Sign orientation(Point_key k, Cell_handle h) const {
    kdel_.set_instantaneous_time();
    int hinf=-1;
    for (int i=0; i< 4; ++i) {
      if (h->vertex(i)->point() == Point_key()) {
	hinf=i;
      }
    }
    CGAL::Sign ret=CGAL::POSITIVE;
    for (int i=0; i< 4; ++i) {
      if (hinf == -1 ||  hinf == i) {
	typename Triangulation::Facet f(h, i);
	typename Traits::Instantaneous_kernel::Orientation_3 o3= triangulation().geom_traits().orientation_3_object();
	
	CGAL::Sign sn= o3((internal::vertex_of_facet(f,0)->point()),
			  (internal::vertex_of_facet(f,1)->point()),
			  (internal::vertex_of_facet(f,2)->point()),
			  (k));
	if (sn ==CGAL::ZERO) {
	  CGAL_KINETIC_LOG(LOG_SOME, "Point " << k << " lies on face ") ;
	  CGAL_KINETIC_LOG_WRITE(LOG_SOME, internal::write_facet( f, LOG_STREAM));
	  CGAL_KINETIC_LOG(LOG_SOME, "\nPoint trajectory is  " << point(k)  << std::endl) ;
	  CGAL_KINETIC_LOG(LOG_SOME, "Triangle 0  " << point(internal::vertex_of_facet(f,0)->point())  << std::endl) ;
	  CGAL_KINETIC_LOG(LOG_SOME, "Triangle 1  " << point(internal::vertex_of_facet(f,1)->point())  << std::endl) ;
	  CGAL_KINETIC_LOG(LOG_SOME, "Triangle 2  " << point(internal::vertex_of_facet(f,2)->point())  << std::endl) ;
	  ret=CGAL::ZERO;
	} else if (sn==CGAL::NEGATIVE) {
	  ret=CGAL::NEGATIVE;
	  return ret;
	}
      }
    }
    return ret;
  }

  void handle_redundant(Point_key k, Cell_handle h, Root_stack s) {
    CGAL_KINETIC_LOG(LOG_LOTS, "Handle redundant " << k << " ") ;
    CGAL_KINETIC_LOG_WRITE(LOG_LOTS, internal::write_cell( h, LOG_STREAM));
    CGAL_KINETIC_LOG(LOG_LOTS, std::endl);
    CGAL_precondition(orientation(k,h) != CGAL::NEGATIVE);
    CGAL_precondition(redundant_points_[k]==Event_key());
    CGAL_assertion_code(bool found=false);
    for( typename RCMap::iterator it = redundant_cells_.begin();  it != redundant_cells_.end(); ++it) {
      CGAL_assertion_code(if(it->second==k) found=true);
      CGAL_assertion(h== it->first || it->second != k);
    }
    CGAL_assertion(found);
    redundant_points_[k]= Event_key();
    if (!kdel_.has_certificates()) return;
 
    Time pst;
    /*if (!ps.empty()) pst = ps.top();
      else pst= std::numeric_limits<Time>::infinity();*/
    if (s.will_fail()) {
      pst=s.failure_time();
    } else {
      pst= kdel_.simulator()->end_time();
    }
    int hinf=-1;
    for (int i=0; i< 4; ++i) {
      if (h->vertex(i)->point() == Point_key()) {
	hinf=i;
      }
    }
    int first=0;
    for (int i=0; i< 4; ++i) {
      if (hinf == -1 || hinf ==i) {
	typename Triangulation::Facet f(h, i);
	// order matters
	Root_stack cs
	  = kdel_.orientation_object()(point(internal::vertex_of_facet(f,0)->point()),
				       point(internal::vertex_of_facet(f,1)->point()),
				       point(internal::vertex_of_facet(f,2)->point()),
				       point(k),
				       kdel_.simulator()->current_time(),
				       kdel_.simulator()->end_time());
	if (cs.will_fail() && cs.failure_time() < pst) {
	  pst= cs.failure_time();
	  s=cs;
	  first= i+1;
	}
      }
    }
    if (pst < kdel_.simulator()->end_time()) {
      s.pop_failure_time();
      if (first==0 ) {
	CGAL_KINETIC_LOG(LOG_LOTS, "Making push certificate for " << k << std::endl);
	redundant_points_[k]= kdel_.simulator()->new_event(pst, Push_event(s, k, h, this));
	CGAL_assertion_code(kdel_.simulator()->audit_event(redundant_points_[k]));
	CGAL_assertion_code(kdel_.simulator()->audit_events());
      } else {
	CGAL_KINETIC_LOG(LOG_LOTS, "Making move certificate for " << k << std::endl);
	redundant_points_[k]= kdel_.simulator()->new_event(pst, Move_event(s, k, h, first-1, this));
	CGAL_assertion_code(kdel_.simulator()->audit_event(redundant_points_[k]));
	CGAL_assertion_code(kdel_.simulator()->audit_events());
      }
    } else {
      redundant_points_[k]= kdel_.simulator()->null_event();
    }

   
  }

  /*
    if (s.will_fail()) {
    Time t= s.failure_time();
    s.pop_failure_time();
    redundant_points_[k]= kdel_.simulator()->new_event(t, Push_event(s, vh, this));
    } else {
    return kdel_.simulator()->null_event();
    }
  */

  void handle_redundant(Point_key k, typename Triangulation::Cell_handle h) {
    int hinf=-1;
    for (int i=0; i< 4; ++i) {
      if (h->vertex(i)->point() == Point_key()) {
	hinf=i;
      }
    }
    Root_stack ps;
    if (hinf ==-1 ) {
      ps
	= kdel_.power_test_object()(point(h->vertex(0)->point()),
				    point(h->vertex(1)->point()),
				    point(h->vertex(2)->point()),
				    point(h->vertex(3)->point()),
				    point(k),
				    kdel_.simulator()->current_time(),
				    kdel_.simulator()->end_time());
    }
    handle_redundant(k,h,ps);
  }


  bool try_handle_redundant(Point_key k, typename Triangulation::Cell_handle h) {
    CGAL_KINETIC_LOG(LOG_LOTS, "Trying handle redundant " << k << " ") ;
    CGAL_KINETIC_LOG_WRITE(LOG_LOTS, internal::write_cell( h, LOG_STREAM));
    CGAL_KINETIC_LOG(LOG_LOTS, std::endl);
    if (orientation(k,h) != CGAL::NEGATIVE) {
      CGAL_precondition(redundant_points_[k]==Event_key());
      redundant_cells_.insert(typename RCMap::value_type(h, k));
      handle_redundant(k,h);
      return true;
    } else {
      return false;
    }
  }

  void handle_vertex(typename Triangulation::Vertex_handle vh, Root_stack &s) {
    CGAL_KINETIC_LOG(LOG_LOTS, "Updating vertex " << vh->point() << std::endl);
    CGAL_precondition(vh->info() == Event_key());
    if (s.will_fail()) {
      Time t= s.failure_time();
      s.pop_failure_time();
      vh->info()= kdel_.simulator()->new_event(t, Pop_event(s, vh, this));
      CGAL_assertion_code(kdel_.simulator()->audit_event(vh->info()));
      CGAL_assertion_code(kdel_.simulator()->audit_events());
      CGAL_assertion_code(kdel_.simulator()->template event<Pop_event>(vh->info()).audit(vh->info()));
    } else {
      vh->info()= kdel_.simulator()->null_event();
    }
  }

  void handle_vertex(typename Triangulation::Vertex_handle vh) {
    CGAL_KINETIC_LOG(LOG_LOTS, "Handling vertex " << vh->point() << std::endl);
    if (vh== triangulation().infinite_vertex()) return;
    CGAL_precondition( internal::has_degree_4(triangulation(), vh));
    CGAL_precondition( vh->info() == Event_key());

    typename Triangulation::Cell_handle ch= vh->cell();
    typename Triangulation::Facet f(ch, ch->index(vh));
    std::vector<typename Triangulation::Vertex_handle> n(4);
   
    for (int i=0; i<3; ++i) {
      n[i]= internal::vertex_of_facet(f,i);
      if (n[i]== triangulation().infinite_vertex()) {
	vh->info() = kdel_.simulator()->null_event();
	return;
      }
    }
    int ind= (f.second+1)%4;
    // some vertex on facet
    n[3] = triangulation().mirror_vertex(ch, ind);
    if (n[3]== triangulation().infinite_vertex()) {
      vh->info() = kdel_.simulator()->null_event();
      return;
    }

    CGAL_KINETIC_LOG(LOG_LOTS, "Making D4 certificate for " << n[0]->point() << n[1]->point()
		     << n[2]->point() << n[3]->point() << " around " << vh->point() << std::endl);

   
    //! The order is switched to invert the predicate since we want it to fail when it goes outside
    Root_stack s
      = kdel_.power_test_object()(point(n[1]->point()),
				  point(n[0]->point()),
				  point(n[2]->point()),
				  point(n[3]->point()),
				  point(vh->point()),
				  kdel_.simulator()->current_time(),
				  kdel_.simulator()->end_time());
    return handle_vertex(vh, s);
  }


  void on_geometry_changed() {
    if (listener_!= NULL) {
      listener_->new_notification(Listener::TRIANGULATION);
    }
    CGAL_KINETIC_LOG(LOG_LOTS, *this);
    audit_structure();
  }

  typename MPT::Data point(Point_key k) const {
    return kdel_.moving_object_table()->at(k);
  }

  // clean vertex events, gather redundant points

  // create vertex events, try to insert redundant points


  void destroy_cell(typename Triangulation::Cell_handle h) {
    CGAL_KINETIC_LOG(LOG_LOTS, "Cleaning cell " << h->vertex(0)->point()
		     << " " << h->vertex(1)->point() << " " << h->vertex(2)->point()
		     << " " << h->vertex(3)->point() << std::endl);
    for (unsigned int i=0; i<4; ++i) {
      if (h->vertex(i)->info() != Event_key()) {
	CGAL_KINETIC_LOG(LOG_LOTS, "Cleaning vertex " << h->vertex(i)->point() << std::endl);
	kdel_.simulator()->delete_event(h->vertex(i)->info());
	h->vertex(i)->info() = Event_key();
      }
    }
    typename RCMap::iterator beg= redundant_cells_.lower_bound(h);
    typename RCMap::iterator end= redundant_cells_.upper_bound(h);
    for (; beg != end; ++beg) {
      unhandled_keys_.push_back(beg->second);
      kdel_.simulator()->delete_event(redundant_points_[beg->second]);
      redundant_points_.erase(beg->second);
    }
    redundant_cells_.erase(redundant_cells_.lower_bound(h),
			   redundant_cells_.upper_bound(h));
  }

  void create_cell(typename Triangulation::Cell_handle h) {
    CGAL_KINETIC_LOG(LOG_LOTS, "Creating cell " << h->vertex(0)->point()
		     << " " << h->vertex(1)->point() << " " << h->vertex(2)->point()
		     << " " << h->vertex(3)->point() << std::endl);
    for (unsigned int i=0; i< 4; ++i){
      if (h->vertex(i)->info() == Event_key() && kdel_.is_degree_4(h->vertex(i))){
	handle_vertex(h->vertex(i));
      }
    }
    for ( int i=0; i< static_cast<int>(unhandled_keys_.size()); ++i) {
      kdel_.set_instantaneous_time(true);
      if (try_handle_redundant(unhandled_keys_[i], h)) {
	unhandled_keys_.erase(unhandled_keys_.begin()+i);
	--i;
      }
    }
  }

  void remove_redundant(typename Triangulation::Cell_handle h, Point_key k) {
    redundant_points_.erase(k);
    typename RCMap::iterator beg = redundant_cells_.lower_bound(h);
    typename RCMap::iterator end = redundant_cells_.upper_bound(h);
    for (; beg != end; ++beg) {
      if (beg->second == k) {
	redundant_cells_.erase(beg);
	return;
      }
    }
    
    for (typename RCMap::iterator cur= redundant_cells_.begin();
	 cur != redundant_cells_.end(); ++cur){
      CGAL_KINETIC_ERROR_WRITE( internal::write_cell( cur->first, LOG_STREAM));
      CGAL_KINETIC_ERROR(": " << cur->second);
      if (cur->second == k) {
	CGAL_assertion_code(Cell_handle ch= cur->first);
	CGAL_assertion(ch==h);
	CGAL_assertion(0);
      }
    }
    CGAL_assertion(0);
  }

 

  KDel kdel_;
  Simulator_listener siml_;
  Moving_point_table_listener motl_;
  Listener *listener_;
  RPMap redundant_points_;
  RCMap redundant_cells_;
  std::vector<Point_key> unhandled_keys_;
  //typename P::Instantaneous_kernel::Orientation_3 po_;
  // typename P::Kinetic_kernel::Weighted_orientation_3 por_;
};

template <class Traits, class Triang, class Visit>
std::ostream &operator<<(std::ostream &out, const Regular_triangulation_3<Traits, Triang, Visit> &rt)
{
  rt.write(out);
  return out;
}


CGAL_KINETIC_END_NAMESPACE

#if defined(BOOST_MSVC)
#  pragma warning(pop)
#endif

#endif