1 /* Copyright (c) 1997-2021
2 Ewgenij Gawrilow, Michael Joswig, and the polymake team
3 Technische Universität Berlin, Germany
4 https://polymake.org
5
6 This program is free software; you can redistribute it and/or modify it
7 under the terms of the GNU General Public License as published by the
8 Free Software Foundation; either version 2, or (at your option) any
9 later version: http://www.gnu.org/licenses/gpl.txt.
10
11 This program is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15 --------------------------------------------------------------------------------
16 */
17
18 #include "polymake/client.h"
19 #include "polymake/Matrix.h"
20 #include "polymake/IncidenceMatrix.h"
21 #include "polymake/Set.h"
22 #include "polymake/topaz/sum_triangulation_tools.h"
23 #include "polymake/PowerSet.h"
24 #include "polymake/hash_set"
25
26 namespace polymake { namespace topaz {
27
28 // this function takes a facet F and glues it around the boundary of
29 // the ball defined by the WEB. It creates a new simplicial complex
30 // which contains all the facets indicated by WEB.
glue_facet(const Set<Int> & F_in,const Array<Int> & F_vertex_indices,const Array<Set<Int>> & facets,const Array<Int> & facets_vertex_indices,const Set<Int> & web,const Int shift,bool shift_facet,std::vector<Set<Int>> & result)31 void glue_facet(const Set<Int>& F_in,
32 const Array<Int>& F_vertex_indices,
33 const Array<Set<Int>>& facets,
34 const Array<Int>& facets_vertex_indices,
35 const Set<Int>& web,
36 const Int shift,
37 bool shift_facet,
38 std::vector<Set<Int>>& result)
39 {
40 // compute the boundary ridges
41 hash_set<Set<Int>> boundary;
42
43 for (const auto& sf : web) {
44 for (auto rit = entire(all_subsets_less_1(facets[sf])); !rit.at_end(); ++rit) {
45 if (boundary.exists(*rit)) {
46 boundary -= *rit;
47 } else {
48 boundary += *rit;
49 }
50 }
51 }
52
53 // take care of index shifting for unused points
54 Array<Int> vertex_indices = facets_vertex_indices;
55
56 // take into account the vertex indices of F
57 Set<Int> F(permuted_inv(F_in, F_vertex_indices));
58
59 // shift the indices of F or boundary facet
60 if (shift_facet) {
61 // TODO: introduce Set::transpose instead of this madness
62 Set<Int> F_shifted(entire(attach_operation(F, same_value(shift), operations::add())));
63 F = F_shifted;
64 } else {
65 // we shift the vertex indices of the boundary facet via the
66 // vertex permutation so that we don't have to do it later
67 for (Int& v : vertex_indices) v+=shift;
68 }
69
70 // glue everything together
71 for (const auto& bf : boundary) {
72 result.push_back(F + permuted_inv(bf, vertex_indices));
73 }
74 }
75
76
77 template <typename Scalar>
sum_triangulation(BigObject p_in,BigObject q_in,const IncidenceMatrix<> webOfStars_in,OptionSet options)78 BigObject sum_triangulation(BigObject p_in,
79 BigObject q_in,
80 const IncidenceMatrix<> webOfStars_in,
81 OptionSet options)
82 {
83 return sum_triangulation_impl<Scalar>(p_in, q_in, webOfStars_in, options);
84 }
85
86 UserFunctionTemplate4perl("# @category Producing a new simplicial complex from others\n"
87 "# Produce a specific sum-triangulation of two given triangulations.\n"
88 "# and a WebOfStars. There are P-sum-triangulations and Q-sum-triangulations."
89 "# If the image of the star of the origin of P is empty then we have a"
90 "# Q-sum-triangulation; otherwise it is a P-sum-triangulation."
91 "# For details see Assarf, Joswig & Pfeifle:"
92 "# Webs of stars or how to triangulate sums of polytopes, to appear"
93 "# @param GeometricSimplicialComplex P first complex"
94 "# @param GeometricSimplicialComplex Q second complex"
95 "# @param IncidenceMatrix WebOfStars Every row corresponds to a full dimensional simplex in P and every column to a full dimensional simplex in Q."
96 "# @option Bool origin_first decides if the origin should be the first point in the resulting complex. Default=0"
97 "# @return GeometricSimplicialComplex",
98 "sum_triangulation<Scalar>(GeometricSimplicialComplex<type_upgrade<Scalar>> GeometricSimplicialComplex<type_upgrade<Scalar>>; IncidenceMatrix=new IncidenceMatrix() { origin_first => 0 })");
99
100 } }
101
102 // Local Variables:
103 // mode:C++
104 // c-basic-offset:3
105 // indent-tabs-mode:nil
106 // End:
107