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 #pragma once
19
20 #include "polymake/Rational.h"
21 #include "polymake/Matrix.h"
22 #include "polymake/TropicalNumber.h"
23
24 namespace polymake{ namespace tropical {
25
26 template <typename Coefficient, typename VType>
27 Vector<Coefficient> thomog_vec(const GenericVector<VType, Coefficient>& affine, Int chart = 0, bool has_leading_coordinate = true)
28 {
29 if (affine.dim() <= 1)
30 return Vector<Coefficient>(affine);
31 if (chart < 0 || chart > (affine.dim()-has_leading_coordinate))
32 throw std::runtime_error("Invalid chart coordinate");
33
34 Vector<Coefficient> proj(affine.dim()+1);
35 proj.slice(~scalar2set(chart+has_leading_coordinate)) = affine;
36 return proj;
37 }
38
39 template <typename Coefficient, typename MType>
40 Matrix<Coefficient> thomog(const GenericMatrix<MType, Coefficient>& affine, Int chart = 0, bool has_leading_coordinate = true)
41 {
42 if (affine.rows() == 0)
43 return Matrix<Coefficient>(0,affine.cols()+1);
44 if (chart < 0 || chart > (affine.cols()-has_leading_coordinate))
45 throw std::runtime_error("Invalid chart coordinate.");
46
47 Matrix<Coefficient> proj(affine.rows(), affine.cols()+1);
48 proj.minor(All,~scalar2set(chart+has_leading_coordinate)) = affine;
49 return proj;
50 }
51
52 template <typename Affine, typename Projective>
tdehomog_elim_col(Affine && affine,const Projective & proj,Int chart,bool has_leading_coordinate)53 void tdehomog_elim_col(Affine&& affine, const Projective& proj, Int chart, bool has_leading_coordinate)
54 {
55 auto elim_col = proj.begin();
56 std::advance(elim_col, chart+has_leading_coordinate);
57 auto it = entire(affine);
58 if (has_leading_coordinate) ++it;
59 for (; !it.at_end(); ++it)
60 *it -= *elim_col;
61 }
62
63 // This does the inverse of thomog, i.e. removes a coordinate by shifting it to 0.
64 template <typename TMatrix, typename Coefficient>
65 Matrix<Coefficient> tdehomog(const GenericMatrix<TMatrix, Coefficient>& proj, Int chart = 0, bool has_leading_coordinate = true)
66 {
67 if (chart < 0 || chart > (proj.cols()-has_leading_coordinate-1))
68 throw std::runtime_error("Invalid chart coordinate");
69
70 Matrix<Coefficient> affine = proj.minor(All, ~scalar2set(chart+has_leading_coordinate));
71 tdehomog_elim_col(cols(affine), cols(proj), chart, has_leading_coordinate);
72 return affine;
73 }
74
75 template <typename TVector, typename Coefficient>
76 Vector<Coefficient> tdehomog_vec(const GenericVector<TVector, Coefficient>& proj, Int chart = 0, bool has_leading_coordinate = true)
77 {
78 if (proj.dim() <= 1)
79 return Vector<Coefficient>();
80 if (chart < 0 || chart > (proj.dim()-has_leading_coordinate-1))
81 throw std::runtime_error("Invalid chart coordinate");
82
83 Vector<Coefficient> affine = proj.slice(~scalar2set(chart+has_leading_coordinate));
84 tdehomog_elim_col(affine, proj.top(), chart, has_leading_coordinate);
85 return affine;
86 }
87
88 /*
89 * @brief: scale the rows of a matrix so that
90 * the first non-null entry of each row is tropical one
91 */
92 template <typename Addition, typename Scalar, typename MatrixTop>
93 Matrix<TropicalNumber<Addition, Scalar>>
normalized_first(const GenericMatrix<MatrixTop,TropicalNumber<Addition,Scalar>> & homogeneous_points)94 normalized_first(const GenericMatrix<MatrixTop, TropicalNumber<Addition, Scalar>>& homogeneous_points)
95 {
96 using TNumber = TropicalNumber<Addition,Scalar>;
97
98 Matrix<TNumber> result(homogeneous_points);
99 for (auto r : rows(result)) {
100 TNumber value;
101 for (auto entry : r) {
102 if (!is_zero(entry)) {
103 value = entry;
104 break;
105 }
106 }
107 if (!is_zero(value)) r /= value;
108 }
109 return result;
110 }
111
112 /*
113 * @brief: scale the rows of a matrix so that
114 * the first non-null entry of each row is tropical one
115 */
116 template <typename Addition, typename Scalar, typename VectorTop>
117 Vector<TropicalNumber<Addition, Scalar>>
normalized_first(const GenericVector<VectorTop,TropicalNumber<Addition,Scalar>> & homogeneous_point)118 normalized_first(const GenericVector<VectorTop, TropicalNumber<Addition, Scalar>>& homogeneous_point)
119 {
120 using TNumber = TropicalNumber<Addition,Scalar>;
121
122 Vector<TNumber> result(homogeneous_point);
123 TNumber value;
124 for (auto entry : result) {
125 if (!is_zero(entry)) {
126 value = entry;
127 break;
128 }
129 }
130 if (!is_zero(value)) result /= value;
131 return result;
132 }
133
134 } }
135
136
137 // Local Variables:
138 // mode:C++
139 // c-basic-offset:3
140 // indent-tabs-mode:nil
141 // End:
142