summaryrefslogtreecommitdiff
path: root/src/Toplex_map/include/gudhi/Toplex_map.h
blob: b6bb5381b46cbf311db6c9d292f5d17ef54f52a1 (plain)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
#ifndef TOPLEX_MAP_H
#define TOPLEX_MAP_H

#include <vector>
#include <unordered_set>
#include <unordered_map>
#include <memory>
#include <limits>

#define vertex_upper_bound std::numeric_limits<Toplex_map::Vertex>::max()

namespace Gudhi {

/** A Toplex_map represents the simplicial complex.
 * A "toplex" is a maximal simplex.
 * \ingroup toplex_map   */
class Toplex_map {

public:

    /** Vertex is the type of vertices.
     * \ingroup toplex_map   */
    typedef std::size_t Vertex;

    /** Simplex is the type of simplices.
     * \ingroup toplex_map   */
    typedef std::unordered_set<Toplex_map::Vertex> Simplex;

    /** The type of the pointers to maximal simplices.
     * \ingroup toplex_map   */
    typedef std::shared_ptr<Toplex_map::Simplex> Simplex_ptr;

    struct Sptr_hash{ std::size_t operator()(const Toplex_map::Simplex_ptr& s) const; };
    struct Sptr_equal{ std::size_t operator()(const Toplex_map::Simplex_ptr& a, const Toplex_map::Simplex_ptr& b) const; };
    /** The type of the sets of Toplex_map::Simplex_ptr.
     * \ingroup toplex_map   */
    typedef std::unordered_set<Toplex_map::Simplex_ptr, Sptr_hash, Sptr_equal> Simplex_ptr_set;

    /** \brief Adds the given simplex to the complex.
     * Nothing happens if the simplex has a coface in the complex.
     * \ingroup toplex_map   */
    template <typename Input_vertex_range>
    void insert_simplex(const Input_vertex_range &vertex_range);

    /** \brief Removes the given simplex and its cofaces from the complex.
     * Its faces are kept inside.
     * \ingroup toplex_map   */
    template <typename Input_vertex_range>
    void remove_simplex(const Input_vertex_range &vertex_range);

    /** Does a simplex belong to the complex ?
     * \ingroup toplex_map   */
    template <typename Input_vertex_range>
    bool membership(const Input_vertex_range &vertex_range) const;

    /** Does a simplex is a toplex ?
     * \ingroup toplex_map   */
    template <typename Input_vertex_range>
    bool maximality(const Input_vertex_range &vertex_range) const;

    /** Gives a set of pointers to the maximal cofaces of a simplex.
     * Gives all the toplices if given the empty simplex.
     * Gives not more than max_number maximal cofaces if max_number is strictly positive.
     * \ingroup toplex_map   */
    template <typename Input_vertex_range>
    Toplex_map::Simplex_ptr_set maximal_cofaces(const Input_vertex_range &vertex_range, const std::size_t max_number = 0) const;

    /** Contracts one edge in the complex.
     * The edge has to verify the link condition if you want to preserve topology.
     * Returns the remaining vertex.
     * \ingroup toplex_map   */
    Toplex_map::Vertex contraction(const Toplex_map::Vertex x, const Toplex_map::Vertex y, bool force=false);

    /** Adds the given simplex to the complex.
     * The simplex must not have neither maximal face nor coface in the complex.
     * \ingroup toplex_map   */
    template <typename Input_vertex_range>
    void insert_independent_simplex(const Input_vertex_range &vertex_range);

    /** \internal Removes a toplex without adding facets after.
     * \ingroup toplex_map   */
    void erase_maximal(const Toplex_map::Simplex_ptr& sptr);

    /** Removes a vertex from any simplex containing it.
     * \ingroup toplex_map   */
    void remove_vertex(const Toplex_map::Vertex x);

    /** \brief Number of maximal simplices.
      * \ingroup toplex_map   */
    std::size_t num_simplices() const;

protected:
    /** \internal Gives an index in order to look for a simplex quickly.
     * \ingroup toplex_map   */
    template <typename Input_vertex_range>
    Toplex_map::Vertex best_index(const Input_vertex_range &vertex_range) const;
    
    /** \internal The map from vertices to toplices
     * \ingroup toplex_map   */
    std::unordered_map<Toplex_map::Vertex, Toplex_map::Simplex_ptr_set> t0;
};

// Pointers are also used as key in the hash sets.
template <typename Input_vertex_range>
Toplex_map::Simplex_ptr get_key(const Input_vertex_range &vertex_range);

// Is the first simplex a face of the second ?
template <typename Input_vertex_range1, typename Input_vertex_range2>
bool included(const Input_vertex_range1 &vertex_range1, const Input_vertex_range2 &vertex_range2);

// All the facets of the given simplex.
template <typename Input_vertex_range>
std::vector<Toplex_map::Simplex> facets(const Input_vertex_range &vertex_range);

template <typename Input_vertex_range>
void Toplex_map::insert_simplex(const Input_vertex_range &vertex_range){
    if(membership(vertex_range)) return;
    bool replace_facets = true;
    for(const Toplex_map::Simplex& facet : facets(vertex_range))
        if(!maximality(facet))
        {
            replace_facets=false;
            break;
        }
    if(replace_facets)
        for(const Toplex_map::Simplex& facet : facets(vertex_range))
            erase_maximal(get_key(facet));
    else
        for(const Toplex_map::Vertex& v : vertex_range)
            if(t0.count(v))  for(const Toplex_map::Simplex_ptr& fptr : Simplex_ptr_set(t0.at(v)))
                //Copy constructor needed because the set is modified
                if(included(*fptr,vertex_range)) erase_maximal(fptr);
    // We erase all the maximal faces of the simplex
    insert_independent_simplex(vertex_range);
}

template <typename Input_vertex_range>
void Toplex_map::remove_simplex(const Input_vertex_range &vertex_range){
    if(vertex_range.begin()==vertex_range.end())
        t0.clear();
    // Removal of the empty simplex means cleaning everything
    else {
        const Toplex_map::Vertex& v = best_index(vertex_range);
        if(t0.count(v)) for(const Toplex_map::Simplex_ptr& sptr : Simplex_ptr_set(t0.at(v)))
            //Copy constructor needed because the set is modified
            if(included(vertex_range, *sptr)){
                erase_maximal(sptr);
                for(const Toplex_map::Simplex& f : facets(vertex_range))
                    if(!membership(f)) insert_independent_simplex(f);
                // We add the facets which are new maximal simplices
            }
    }
}

template <typename Input_vertex_range>
bool Toplex_map::membership(const Input_vertex_range &vertex_range) const{
    if(t0.size()==0) return false;
    const Toplex_map::Vertex& v = best_index(vertex_range);
    if(!t0.count(v))  return false;
    if(maximality(vertex_range)) return true;
    for(const Toplex_map::Simplex_ptr& sptr : t0.at(v))
        if(included(vertex_range, *sptr))
            return true;
    return false;
}

template <typename Input_vertex_range>
bool Toplex_map::maximality(const Input_vertex_range &vertex_range) const{
    const Toplex_map::Vertex& v =  best_index(vertex_range);
    if(!t0.count(v)) return false;
    return t0.at(v).count(get_key(vertex_range));
}

template <typename Input_vertex_range>
Toplex_map::Simplex_ptr_set Toplex_map::maximal_cofaces(const Input_vertex_range &vertex_range, const std::size_t max_number) const{
    Simplex_ptr_set cofaces;
    if(maximality(vertex_range))
        cofaces.emplace(get_key(vertex_range));
    else if(vertex_range.begin()==vertex_range.end())
        for(const auto& kv : t0)
            for(const Toplex_map::Simplex_ptr& sptr : kv.second){
                //kv.second is a Simplex_ptr_set
                cofaces.emplace(sptr);
                if(cofaces.size()==max_number)
                    return cofaces;
            }
    else {
        const Toplex_map::Vertex& v = best_index(vertex_range);
        if(t0.count(v)) for(const Toplex_map::Simplex_ptr& sptr : t0.at(v))
            if(included(vertex_range, *sptr)){
                cofaces.emplace(sptr);
                if(cofaces.size()==max_number)
                    return cofaces;
            }
    }
    return cofaces;
}

Toplex_map::Vertex Toplex_map::contraction(const Toplex_map::Vertex x, const Toplex_map::Vertex y, bool force){
    if(!t0.count(x)) return y;
    if(!t0.count(y)) return x;
    int k, d;
    if(force || (t0.at(x).size() > t0.at(y).size()))
        k=x, d=y;
    else
        k=y, d=x;
    for(const Toplex_map::Simplex_ptr& sptr : Simplex_ptr_set(t0.at(d))){
        //Copy constructor needed because the set is modified
        Simplex sigma(*sptr);
        Simplex s; s.insert(2);
        erase_maximal(sptr);
        sigma.erase(d);
        sigma.insert(k);
        insert_simplex(sigma);
    }
    return k;
}

template <typename Input_vertex_range>
void Toplex_map::insert_independent_simplex(const Input_vertex_range &vertex_range){
    for(const Toplex_map::Vertex& v : vertex_range){
        if(!t0.count(v)) t0.emplace(v, Simplex_ptr_set());
        auto k = get_key(vertex_range);
        t0.at(v).emplace(k);
    }
}

void Toplex_map::remove_vertex(const Toplex_map::Vertex x){
    for(const Toplex_map::Simplex_ptr& sptr : Simplex_ptr_set(t0.at(x))){
        Simplex sigma(*sptr);
        erase_maximal(sptr);
        sigma.erase(x);
        insert_simplex(sigma);
    }
}

std::size_t Toplex_map::num_simplices() const{
    return maximal_cofaces(Simplex()).size();
}

inline void Toplex_map::erase_maximal(const Toplex_map::Simplex_ptr& sptr){
    Simplex sigma(*sptr);
    if (sptr->size()==0)
        sigma.insert(vertex_upper_bound);
    for(const Toplex_map::Vertex& v : sigma){
        t0.at(v).erase(sptr);
        if(t0.at(v).size()==0) t0.erase(v);
    }
}

template <typename Input_vertex_range>
Toplex_map::Vertex Toplex_map::best_index(const Input_vertex_range &vertex_range) const{
    std::size_t min = std::numeric_limits<size_t>::max();
    Vertex arg_min = vertex_upper_bound;
    for(const Toplex_map::Vertex& v : vertex_range)
        if(!t0.count(v)) return v;
        else if(t0.at(v).size() < min)
            min = t0.at(v).size(), arg_min = v;
    return arg_min;
}

std::size_t Toplex_map::Sptr_equal::operator()(const Toplex_map::Simplex_ptr& s1, const Toplex_map::Simplex_ptr& s2) const {
    if (s1->size() != s2->size()) return false;
    return included(*s1,*s2);
    // inclusion tests equality for same size simplices
}

std::size_t Toplex_map::Sptr_hash::operator()(const Toplex_map::Simplex_ptr& s) const {
    std::hash<double> h_f;
    //double hash works better than int hash
    size_t h = 0;
    for(const Toplex_map::Vertex& v : *s)
        h += h_f(static_cast<double>(v));
    return h;
}

template <typename Input_vertex_range>
Toplex_map::Simplex_ptr get_key(const Input_vertex_range &vertex_range){
    Toplex_map::Simplex s(vertex_range.begin(), vertex_range.end());
    return std::make_shared<Toplex_map::Simplex>(s);
}

template <typename Input_vertex_range1, typename Input_vertex_range2>
bool included(const Input_vertex_range1 &vertex_range1, const Input_vertex_range2 &vertex_range2){
    Toplex_map::Simplex s2(vertex_range2.begin(), vertex_range2.end());
    for(const Toplex_map::Vertex& v : vertex_range1)
        if(!s2.count(v)) return false;
    return true;
}

template <typename Input_vertex_range>
std::vector<Toplex_map::Simplex> facets(const Input_vertex_range &vertex_range){
    std::vector<Toplex_map::Simplex> facets;
    Toplex_map::Simplex f(vertex_range.begin(), vertex_range.end());
    for(const Toplex_map::Vertex& v : vertex_range){
        f.erase(v);
        facets.emplace_back(f);
        f.insert(v);
    }
    return facets;
}

} //namespace Gudhi

#endif /* TOPLEX_MAP_H */