123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483 |
- //=======================================================================
- // Copyright 1997, 1998, 1999, 2000 University of Notre Dame.
- // Authors: Andrew Lumsdaine, Lie-Quan Lee, Jeremy G. Siek
- //
- // Distributed under the Boost Software License, Version 1.0. (See
- // accompanying file LICENSE_1_0.txt or copy at
- // http://www.boost.org/LICENSE_1_0.txt)
- //=======================================================================
- #ifndef BOOST_SELF_AVOIDING_WALK_HPP
- #define BOOST_SELF_AVOIDING_WALK_HPP
- /*
- This file defines necessary components for SAW.
- mesh language: (defined by myself to clearify what is what)
- A triangle in mesh is called an triangle.
- An edge in mesh is called an line.
- A vertex in mesh is called a point.
- A triangular mesh corresponds to a graph in which a vertex is a
- triangle and an edge(u, v) stands for triangle u and triangle v
- share an line.
- After this point, a vertex always refers to vertex in graph,
- therefore it is a traingle in mesh.
- */
- #include <utility>
- #include <boost/config.hpp>
- #include <boost/graph/graph_traits.hpp>
- #include <boost/property_map/property_map.hpp>
- #define SAW_SENTINAL -1
- namespace boost
- {
- template < class T1, class T2, class T3 > struct triple
- {
- T1 first;
- T2 second;
- T3 third;
- triple(const T1& a, const T2& b, const T3& c)
- : first(a), second(b), third(c)
- {
- }
- triple() : first(SAW_SENTINAL), second(SAW_SENTINAL), third(SAW_SENTINAL) {}
- };
- typedef triple< int, int, int > Triple;
- /* Define a vertex property which has a triangle inside. Triangle is
- represented by a triple. */
- struct triangle_tag
- {
- enum
- {
- num = 100
- };
- };
- typedef property< triangle_tag, Triple > triangle_property;
- /* Define an edge property with a line. A line is represented by a
- pair. This is not required for SAW though.
- */
- struct line_tag
- {
- enum
- {
- num = 101
- };
- };
- template < class T >
- struct line_property : public property< line_tag, std::pair< T, T > >
- {
- };
- /*Precondition: Points in a Triangle are in order */
- template < class Triangle, class Line >
- inline void get_sharing(const Triangle& a, const Triangle& b, Line& l)
- {
- l.first = SAW_SENTINAL;
- l.second = SAW_SENTINAL;
- if (a.first == b.first)
- {
- l.first = a.first;
- if (a.second == b.second || a.second == b.third)
- l.second = a.second;
- else if (a.third == b.second || a.third == b.third)
- l.second = a.third;
- }
- else if (a.first == b.second)
- {
- l.first = a.first;
- if (a.second == b.third)
- l.second = a.second;
- else if (a.third == b.third)
- l.second = a.third;
- }
- else if (a.first == b.third)
- {
- l.first = a.first;
- }
- else if (a.second == b.first)
- {
- l.first = a.second;
- if (a.third == b.second || a.third == b.third)
- l.second = a.third;
- }
- else if (a.second == b.second)
- {
- l.first = a.second;
- if (a.third == b.third)
- l.second = a.third;
- }
- else if (a.second == b.third)
- {
- l.first = a.second;
- }
- else if (a.third == b.first || a.third == b.second || a.third == b.third)
- l.first = a.third;
- /*Make it in order*/
- if (l.first > l.second)
- {
- typename Line::first_type i = l.first;
- l.first = l.second;
- l.second = i;
- }
- }
- template < class TriangleDecorator, class Vertex, class Line >
- struct get_vertex_sharing
- {
- typedef std::pair< Vertex, Line > Pair;
- get_vertex_sharing(const TriangleDecorator& _td) : td(_td) {}
- inline Line operator()(const Vertex& u, const Vertex& v) const
- {
- Line l;
- get_sharing(td[u], td[v], l);
- return l;
- }
- inline Line operator()(const Pair& u, const Vertex& v) const
- {
- Line l;
- get_sharing(td[u.first], td[v], l);
- return l;
- }
- inline Line operator()(const Pair& u, const Pair& v) const
- {
- Line l;
- get_sharing(td[u.first], td[v.first], l);
- return l;
- }
- TriangleDecorator td;
- };
- /* HList has to be a handle of data holder so that pass-by-value is
- * in right logic.
- *
- * The element of HList is a pair of vertex and line. (remember a
- * line is a pair of two ints.). That indicates the walk w from
- * current vertex is across line. (If the first of line is -1, it is
- * a point though.
- */
- template < class TriangleDecorator, class HList, class IteratorD >
- class SAW_visitor : public bfs_visitor<>, public dfs_visitor<>
- {
- typedef typename boost::property_traits< IteratorD >::value_type iter;
- /*use boost shared_ptr*/
- typedef typename HList::element_type::value_type::second_type Line;
- public:
- typedef tree_edge_tag category;
- inline SAW_visitor(TriangleDecorator _td, HList _hlist, IteratorD ia)
- : td(_td), hlist(_hlist), iter_d(ia)
- {
- }
- template < class Vertex, class Graph >
- inline void start_vertex(Vertex v, Graph&)
- {
- Line l1;
- l1.first = SAW_SENTINAL;
- l1.second = SAW_SENTINAL;
- hlist->push_front(std::make_pair(v, l1));
- iter_d[v] = hlist->begin();
- }
- /*Several symbols:
- w(i): i-th triangle in walk w
- w(i) |- w(i+1): w enter w(i+1) from w(i) over a line
- w(i) ~> w(i+1): w enter w(i+1) from w(i) over a point
- w(i) -> w(i+1): w enter w(i+1) from w(i)
- w(i) ^ w(i+1): the line or point w go over from w(i) to w(i+1)
- */
- template < class Edge, class Graph > bool tree_edge(Edge e, Graph& G)
- {
- using std::make_pair;
- typedef typename boost::graph_traits< Graph >::vertex_descriptor Vertex;
- Vertex tau = target(e, G);
- Vertex i = source(e, G);
- get_vertex_sharing< TriangleDecorator, Vertex, Line > get_sharing_line(
- td);
- Line tau_i = get_sharing_line(tau, i);
- iter w_end = hlist->end();
- iter w_i = iter_d[i];
- iter w_i_m_1 = w_i;
- iter w_i_p_1 = w_i;
- /*----------------------------------------------------------
- * true false
- *==========================================================
- *a w(i-1) |- w(i) w(i-1) ~> w(i) or w(i-1) is null
- *----------------------------------------------------------
- *b w(i) |- w(i+1) w(i) ~> w(i+1) or no w(i+1) yet
- *----------------------------------------------------------
- */
- bool a = false, b = false;
- --w_i_m_1;
- ++w_i_p_1;
- b = (w_i->second.first != SAW_SENTINAL);
- if (w_i_m_1 != w_end)
- {
- a = (w_i_m_1->second.first != SAW_SENTINAL);
- }
- if (a)
- {
- if (b)
- {
- /*Case 1:
- w(i-1) |- w(i) |- w(i+1)
- */
- Line l1 = get_sharing_line(*w_i_m_1, tau);
- iter w_i_m_2 = w_i_m_1;
- --w_i_m_2;
- bool c = true;
- if (w_i_m_2 != w_end)
- {
- c = w_i_m_2->second != l1;
- }
- if (c)
- { /* w(i-1) ^ tau != w(i-2) ^ w(i-1) */
- /*extension: w(i-1) -> tau |- w(i) */
- w_i_m_1->second = l1;
- /*insert(pos, const T&) is to insert before pos*/
- iter_d[tau] = hlist->insert(w_i, make_pair(tau, tau_i));
- }
- else
- { /* w(i-1) ^ tau == w(i-2) ^ w(i-1) */
- /*must be w(i-2) ~> w(i-1) */
- bool d = true;
- // need to handle the case when w_i_p_1 is null
- Line l3 = get_sharing_line(*w_i_p_1, tau);
- if (w_i_p_1 != w_end)
- d = w_i_p_1->second != l3;
- if (d)
- { /* w(i+1) ^ tau != w(i+1) ^ w(i+2) */
- /*extension: w(i) |- tau -> w(i+1) */
- w_i->second = tau_i;
- iter_d[tau]
- = hlist->insert(w_i_p_1, make_pair(tau, l3));
- }
- else
- { /* w(i+1) ^ tau == w(i+1) ^ w(i+2) */
- /*must be w(1+1) ~> w(i+2) */
- Line l5 = get_sharing_line(*w_i_m_1, *w_i_p_1);
- if (l5 != w_i_p_1->second)
- { /* w(i-1) ^ w(i+1) != w(i+1) ^ w(i+2) */
- /*extension: w(i-2) -> tau |- w(i) |- w(i-1) ->
- * w(i+1) */
- w_i_m_2->second = get_sharing_line(*w_i_m_2, tau);
- iter_d[tau]
- = hlist->insert(w_i, make_pair(tau, tau_i));
- w_i->second = w_i_m_1->second;
- w_i_m_1->second = l5;
- iter_d[w_i_m_1->first]
- = hlist->insert(w_i_p_1, *w_i_m_1);
- hlist->erase(w_i_m_1);
- }
- else
- {
- /*mesh is tetrahedral*/
- // dont know what that means.
- ;
- }
- }
- }
- }
- else
- {
- /*Case 2:
- w(i-1) |- w(i) ~> w(1+1)
- */
- if (w_i->second.second == tau_i.first
- || w_i->second.second == tau_i.second)
- { /*w(i) ^ w(i+1) < w(i) ^ tau*/
- /*extension: w(i) |- tau -> w(i+1) */
- w_i->second = tau_i;
- Line l1 = get_sharing_line(*w_i_p_1, tau);
- iter_d[tau] = hlist->insert(w_i_p_1, make_pair(tau, l1));
- }
- else
- { /*w(i) ^ w(i+1) !< w(i) ^ tau*/
- Line l1 = get_sharing_line(*w_i_m_1, tau);
- bool c = true;
- iter w_i_m_2 = w_i_m_1;
- --w_i_m_2;
- if (w_i_m_2 != w_end)
- c = l1 != w_i_m_2->second;
- if (c)
- { /*w(i-1) ^ tau != w(i-2) ^ w(i-1)*/
- /*extension: w(i-1) -> tau |- w(i)*/
- w_i_m_1->second = l1;
- iter_d[tau] = hlist->insert(w_i, make_pair(tau, tau_i));
- }
- else
- { /*w(i-1) ^ tau == w(i-2) ^ w(i-1)*/
- /*must be w(i-2)~>w(i-1)*/
- /*extension: w(i-2) -> tau |- w(i) |- w(i-1) -> w(i+1)*/
- w_i_m_2->second = get_sharing_line(*w_i_m_2, tau);
- iter_d[tau] = hlist->insert(w_i, make_pair(tau, tau_i));
- w_i->second = w_i_m_1->second;
- w_i_m_1->second = get_sharing_line(*w_i_m_1, *w_i_p_1);
- iter_d[w_i_m_1->first]
- = hlist->insert(w_i_p_1, *w_i_m_1);
- hlist->erase(w_i_m_1);
- }
- }
- }
- }
- else
- {
- if (b)
- {
- /*Case 3:
- w(i-1) ~> w(i) |- w(i+1)
- */
- bool c = false;
- if (w_i_m_1 != w_end)
- c = (w_i_m_1->second.second == tau_i.first)
- || (w_i_m_1->second.second == tau_i.second);
- if (c)
- { /*w(i-1) ^ w(i) < w(i) ^ tau*/
- /* extension: w(i-1) -> tau |- w(i) */
- if (w_i_m_1 != w_end)
- w_i_m_1->second = get_sharing_line(*w_i_m_1, tau);
- iter_d[tau] = hlist->insert(w_i, make_pair(tau, tau_i));
- }
- else
- {
- bool d = true;
- Line l1;
- l1.first = SAW_SENTINAL;
- l1.second = SAW_SENTINAL;
- if (w_i_p_1 != w_end)
- {
- l1 = get_sharing_line(*w_i_p_1, tau);
- d = l1 != w_i_p_1->second;
- }
- if (d)
- { /*w(i+1) ^ tau != w(i+1) ^ w(i+2)*/
- /*extension: w(i) |- tau -> w(i+1) */
- w_i->second = tau_i;
- iter_d[tau]
- = hlist->insert(w_i_p_1, make_pair(tau, l1));
- }
- else
- {
- /*must be w(i+1) ~> w(i+2)*/
- /*extension: w(i-1) -> w(i+1) |- w(i) |- tau -> w(i+2)
- */
- iter w_i_p_2 = w_i_p_1;
- ++w_i_p_2;
- w_i_p_1->second = w_i->second;
- iter_d[i] = hlist->insert(w_i_p_2, make_pair(i, tau_i));
- hlist->erase(w_i);
- Line l2 = get_sharing_line(*w_i_p_2, tau);
- iter_d[tau]
- = hlist->insert(w_i_p_2, make_pair(tau, l2));
- }
- }
- }
- else
- {
- /*Case 4:
- w(i-1) ~> w(i) ~> w(i+1)
- */
- bool c = false;
- if (w_i_m_1 != w_end)
- {
- c = (w_i_m_1->second.second == tau_i.first)
- || (w_i_m_1->second.second == tau_i.second);
- }
- if (c)
- { /*w(i-1) ^ w(i) < w(i) ^ tau */
- /*extension: w(i-1) -> tau |- w(i) */
- if (w_i_m_1 != w_end)
- w_i_m_1->second = get_sharing_line(*w_i_m_1, tau);
- iter_d[tau] = hlist->insert(w_i, make_pair(tau, tau_i));
- }
- else
- {
- /*extension: w(i) |- tau -> w(i+1) */
- w_i->second = tau_i;
- Line l1;
- l1.first = SAW_SENTINAL;
- l1.second = SAW_SENTINAL;
- if (w_i_p_1 != w_end)
- l1 = get_sharing_line(*w_i_p_1, tau);
- iter_d[tau] = hlist->insert(w_i_p_1, make_pair(tau, l1));
- }
- }
- }
- return true;
- }
- protected:
- TriangleDecorator td; /*a decorator for vertex*/
- HList hlist;
- /*This must be a handle of list to record the SAW
- The element type of the list is pair<Vertex, Line>
- */
- IteratorD iter_d;
- /*Problem statement: Need a fast access to w for triangle i.
- *Possible solution: mantain an array to record.
- iter_d[i] will return an iterator
- which points to w(i), where i is a vertex
- representing triangle i.
- */
- };
- template < class Triangle, class HList, class Iterator >
- inline SAW_visitor< Triangle, HList, Iterator > visit_SAW(
- Triangle t, HList hl, Iterator i)
- {
- return SAW_visitor< Triangle, HList, Iterator >(t, hl, i);
- }
- template < class Tri, class HList, class Iter >
- inline SAW_visitor< random_access_iterator_property_map< Tri*, Tri, Tri& >,
- HList, random_access_iterator_property_map< Iter*, Iter, Iter& > >
- visit_SAW_ptr(Tri* t, HList hl, Iter* i)
- {
- typedef random_access_iterator_property_map< Tri*, Tri, Tri& > TriD;
- typedef random_access_iterator_property_map< Iter*, Iter, Iter& > IterD;
- return SAW_visitor< TriD, HList, IterD >(t, hl, i);
- }
- // should also have combo's of pointers, and also const :(
- }
- #endif /*BOOST_SAW_H*/
|