graph_subgraph_isomorphism.hh 11.1 KB
Newer Older
1
2
// graph-tool -- a general graph modification and manipulation thingy
//
Tiago Peixoto's avatar
Tiago Peixoto committed
3
// Copyright (C) 2006-2014 Tiago de Paula Peixoto <tiago@skewed.de>
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
//
// This program is free software; you can redistribute it and/or
// modify it under the terms of the GNU General Public License
// as published by the Free Software Foundation; either version 3
// of the License, or (at your option) any later version.
//
// This program is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
// GNU General Public License for more details.
//
// You should have received a copy of the GNU General Public License
// along with this program. If not, see <http://www.gnu.org/licenses/>.

#ifndef SUBGRAPH_ISOMORPHISM_HPP
#define SUBGRAPH_ISOMORPHISM_HPP

#include <boost/graph/graph_traits.hpp>
#include <utility>
Tiago Peixoto's avatar
Tiago Peixoto committed
23
#include <unordered_set>
24
25
26
27
28

namespace boost
{
using namespace std;

29

30
31
32
namespace detail {

//sparse matrix
Tiago Peixoto's avatar
Tiago Peixoto committed
33
typedef vector<std::unordered_set<size_t> > matrix_t;
34
35
36
37
38
39

struct check_adjacency
{
    template <class Graph>
    typename graph_traits<Graph>::vertex_descriptor
    get_vertex(const typename graph_traits<Graph>::edge_descriptor& e, Graph& g,
Tiago Peixoto's avatar
Tiago Peixoto committed
40
               std::true_type out_edges)
41
42
43
44
45
46
47
    {
        return target(e, g);
    }

    template <class Graph>
    typename graph_traits<Graph>::vertex_descriptor
    get_vertex(const typename graph_traits<Graph>::edge_descriptor& e, Graph& g,
Tiago Peixoto's avatar
Tiago Peixoto committed
48
               std::false_type out_edges)
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
    {
        return source(e, g);
    }

    template <class Graph, class IsOut>
    struct get_edge_iterator
    {
        typedef typename graph_traits<Graph>::out_edge_iterator type;

        static pair<type,type>
        edges(typename graph_traits<Graph>::vertex_descriptor v, Graph& g)
        {
            return out_edges(v, g);
        }
    };

    template <class Graph>
Tiago Peixoto's avatar
Tiago Peixoto committed
66
    struct get_edge_iterator<Graph, std::false_type>
67
68
69
70
71
72
73
74
75
76
77
78
79
80
    {
        typedef typename graph_traits<Graph>::in_edge_iterator type;

        static pair<type,type>
        edges(typename graph_traits<Graph>::vertex_descriptor v, Graph& g)
        {
            return in_edges(v, g);
        }
    };

    template <class Graph1, class Graph2, class EdgeLabelling>
    bool operator()(typename graph_traits<Graph1>::vertex_descriptor k,
                    typename graph_traits<Graph2>::vertex_descriptor l,
                    matrix_t& M, EdgeLabelling& edge_labelling,
81
                    Graph1& g1, Graph2& g2, vector<size_t>& vindex,
Tiago Peixoto's avatar
Tiago Peixoto committed
82
                    std::true_type directed)
83
    {
84
        return do_check(k, l, M, edge_labelling, g1, g2, vindex,
Tiago Peixoto's avatar
Tiago Peixoto committed
85
86
                        std::true_type()) &&
            do_check(k, l, M, edge_labelling, g1, g2, vindex, std::false_type());
87
88
89
90
91
    }

    template <class Graph1, class Graph2, class EdgeLabelling>
    bool operator()(typename graph_traits<Graph1>::vertex_descriptor k,
                    typename graph_traits<Graph2>::vertex_descriptor l,
92
                    matrix_t& M, EdgeLabelling& edge_labelling, Graph1& g1,
Tiago Peixoto's avatar
Tiago Peixoto committed
93
                    Graph2& g2, vector<size_t>& vindex, std::false_type directed)
94
    {
Tiago Peixoto's avatar
Tiago Peixoto committed
95
        return do_check(k, l, M, edge_labelling, g1, g2, vindex, std::true_type());
96
97
98
99
100
101
    }

    template <class Graph1, class Graph2, class EdgeLabelling, class IsOut>
    bool do_check(typename graph_traits<Graph1>::vertex_descriptor k,
                  typename graph_traits<Graph2>::vertex_descriptor l,
                  matrix_t& M, EdgeLabelling& edge_labelling, Graph1& g1,
102
                  Graph2& g2, vector<size_t>& vindex, IsOut)
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
    {
        bool valid = true;
        typename get_edge_iterator<Graph1, IsOut>::type e1, e1_end;
        for (tie(e1, e1_end) =
                 get_edge_iterator<Graph1, IsOut>::edges(k, g1);
             e1 != e1_end; ++e1)
        {
            typename graph_traits<Graph1>::vertex_descriptor v1 =
                get_vertex(*e1, g1, IsOut());

            bool is_adjacent = false;
            typename get_edge_iterator<Graph2, IsOut>::type e2, e2_end;
            for (tie(e2, e2_end) =
                     get_edge_iterator<Graph2, IsOut>::edges(l, g2);
                 e2 != e2_end; ++e2)
            {
                typename graph_traits<Graph2>::vertex_descriptor v2 =
                    get_vertex(*e2, g2, IsOut());

122
123
                if (M[v1].find(vindex[v2]) != M[v1].end() &&
                    edge_labelling(*e1, *e2))
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
                {
                    is_adjacent = true;
                    break;
                }
            }

            if (!is_adjacent)
            {
                valid = false;
                break;
            }
        }
        return valid;
    }
};


template <class Graph1, class Graph2, class EdgeLabelling>
142
bool refine_check(const Graph1& sub, const Graph2& g, matrix_t& M, size_t count,
Tiago Peixoto's avatar
Tiago Peixoto committed
143
                  std::unordered_set<size_t>& already_mapped,
144
145
                  EdgeLabelling edge_labelling, vector<size_t>& vlist,
                  vector<size_t>& vindex)
146
{
147
    matrix_t M_temp(num_vertices(sub));
148

149
    int k = 0, N = num_vertices(sub);
150
    #pragma omp parallel for default(shared) private(k) schedule(runtime) if (N > 100)
151
152
153
154
155
156
157
158
    for (k = 0; k < int(count); ++k)
        M_temp[k] = M[k];

    size_t n_mod = 1;
    while (n_mod > 0)
    {
        n_mod = 0;
        bool abort = false;
159
        #pragma omp parallel for default(shared) private(k) schedule(runtime) if (N > 100) \
160
161
162
163
164
            reduction(+:n_mod)
        for (k = count; k < N; ++k)
        {
            if (abort)
                continue;
165
            if (vertex(k, sub) == graph_traits<Graph1>::null_vertex())
166
                continue;
Tiago Peixoto's avatar
Tiago Peixoto committed
167
            std::unordered_set<size_t> m_new;
168
169
170
171
172
            for (typeof(M[k].begin()) li = M[k].begin(); li != M[k].end(); ++li)
            {
                size_t l = *li;
                if (already_mapped.find(l) != already_mapped.end())
                    continue;
173
174
175
176
                bool valid = check_adjacency()(vertex(k, sub),
                                               vertex(vlist[l], g), M,
                                               edge_labelling, sub, g, vindex,
                                               typename graph_tool::is_directed::apply<Graph1>::type());
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
                if (valid)
                    m_new.insert(l);
            }
            if (m_new.empty())
            {
                abort = true;
                continue;
            }
            M_temp[k].swap(m_new);
            if (M_temp[k].size() < M[k].size())
                n_mod++;
        }
        if (abort)
            return false;
        M.swap(M_temp);
    }
    return true;
}


template <class Graph1, class Graph2, class EdgeLabelling, class Mapping>
198
void find_mappings(const Graph1& sub, const Graph2& g, matrix_t& M0,
199
                   vector<Mapping>& FF, EdgeLabelling edge_labelling,
200
                   vector<size_t>& vlist, vector<size_t>& vindex, size_t max_n)
201
202
{
    size_t i = 0;
203
204
    for (i = 0; i < num_vertices(sub); ++i)
        if (vertex(i, sub) != graph_traits<Graph1>::null_vertex())
205
            break;
206
207
208
    int last_i = num_vertices(sub) - 1;
    for (; last_i >= 0; --last_i)
        if (vertex(last_i, sub) != graph_traits<Graph1>::null_vertex())
209
210
211
            break;

    Mapping F;
212
    // [current M] [current sub vertex] [current mapping vertex]
Tiago Peixoto's avatar
Tiago Peixoto committed
213
214
    typedef std::tuple<matrix_t, size_t,
                       typename matrix_t::value_type::const_iterator> state_t;
215
    std::list<state_t> Mstack;
Tiago Peixoto's avatar
Tiago Peixoto committed
216
217
218
    Mstack.push_back(std::make_tuple(M0, i, M0[i].begin()));
    std::get<2>(Mstack.back()) = std::get<0>(Mstack.back())[i].begin();
    std::unordered_set<size_t> already_mapped;
219
220

    // perform depth-first search of combination space
221
    while (!Mstack.empty() && (max_n == 0 || FF.size() < max_n))
222
    {
223
        state_t& state = Mstack.back();
Tiago Peixoto's avatar
Tiago Peixoto committed
224
225
226
        const matrix_t& M = std::get<0>(state);
        size_t& i = std::get<1>(state);
        typename matrix_t::value_type::const_iterator& mi = std::get<2>(state);
227

228
        if (mi == M[i].end()) // no more candidate mappings available
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
        {
            Mstack.pop_back();
            if (!F.empty())
            {
                already_mapped.erase(F.back().second);
                F.pop_back();
            }
            continue;
        }

        matrix_t M_prime(M);
        M_prime[i].clear();
        M_prime[i].insert(*mi);
        already_mapped.insert(*mi);
        size_t c_mi = *mi;

        // move locally to next child
        ++mi;

        size_t ni = i + 1;
249
250
        for (; ni < num_vertices(sub); ++ni)
            if (vertex(ni, sub) != graph_traits<Graph1>::null_vertex())
251
252
253
                break;

        // refine search tree
254
255
        if (refine_check(sub, g, M_prime, ni, already_mapped, edge_labelling,
                         vlist, vindex))
256
257
258
259
260
261
262
        {
            // store the current mapping so far
            F.push_back(std::make_pair(i, c_mi));

            if (ni < size_t(last_i))
            {
                // proceed with search at a higher depth
Tiago Peixoto's avatar
Tiago Peixoto committed
263
264
                Mstack.push_back(std::make_tuple(M_prime, ni, M_prime[ni].begin()));
                std::get<2>(Mstack.back()) = std::get<0>(Mstack.back())[ni].begin();
265
266
267
268
269
270
271
272
            }
            else
            {
                // maximum depth reached: visit all end leafs
                for (typeof(M_prime[ni].begin()) iter = M_prime[ni].begin();
                     iter != M_prime[ni].end(); ++iter)
                {
                    F.push_back(std::make_pair(ni, *iter));
273
274
                    if (max_n == 0 || FF.size() < max_n)
                        FF.push_back(F);
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
                    F.pop_back();
                }

                // we are done which this tree node
                F.pop_back();
                already_mapped.erase(c_mi);
            }
        }
        else
        {
            already_mapped.erase(c_mi);
        }
    }
}

}  // namespace detail


template <class Graph1, class Graph2, class VertexLabelling,
          class EdgeLabelling, class Mapping>
295
void subgraph_isomorphism(const Graph1& sub, const Graph2& g,
296
                          VertexLabelling vertex_labelling,
297
                          EdgeLabelling edge_labelling, vector<Mapping>& F,
298
                          vector<size_t>& vlist, size_t max_n)
299
{
300
301
302
303
304
    // initial mapping candidates
    detail::matrix_t M0(num_vertices(sub));
    vector<size_t> vindex(num_vertices(g));
    for (size_t j = 0; j < num_vertices(g); ++j)
        vindex[vlist[j]] = j;
305
306

    bool abort = false;
307
    int i, N = num_vertices(sub);
308
    #pragma omp parallel for default(shared) private(i) schedule(runtime) if (N > 100)
309
310
    for (i = 0; i < N; ++i)
    {
311
        if (vertex(i, sub) == graph_traits<Graph1>::null_vertex() || abort)
312
313
            continue;

314
        for (size_t j = 0; j < num_vertices(g); ++j)
315
        {
316
            if (vertex(vlist[j], g) == graph_traits<Graph1>::null_vertex())
317
                continue;
318
            if (vertex_labelling(vertex(i, sub), vertex(vlist[j], g)))
319
320
321
322
323
324
325
                M0[i].insert(j);
        }
        if (M0[i].empty())
            abort = true;
    }
    if (abort)
        return;
326
    detail::find_mappings(sub, g, M0, F, edge_labelling, vlist, vindex, max_n);
327
328
329
330
331
332
}

}  // namespace boost


#endif // SUBGRAPH_ISOMORPHISM_HPP