1 // Copyright (c) 2018 ERGO-Code. See license.txt for license.
2
3 #include "utils.h"
4 #include <algorithm>
5 #include <cassert>
6 #include <cmath>
7 #include <utility>
8
9 namespace ipx {
10
AllFinite(const Vector & x)11 bool AllFinite(const Vector& x) {
12 for (double xi : x)
13 if (!std::isfinite(xi))
14 return false;
15 return true;
16 }
17
Onenorm(const Vector & x)18 double Onenorm(const Vector& x) {
19 double norm = 0.0;
20 for (double xi : x)
21 norm += std::abs(xi);
22 return norm;
23 }
24
Twonorm(const Vector & x)25 double Twonorm(const Vector& x) {
26 double norm = 0.0;
27 for (double xi : x)
28 norm += xi * xi;
29 return std::sqrt(norm);
30 }
31
Infnorm(const Vector & x)32 double Infnorm(const Vector& x) {
33 double norm = 0.0;
34 for (double xi : x)
35 norm = std::max(norm, std::abs(xi));
36 return norm;
37 }
38
Dot(const Vector & x,const Vector & y)39 double Dot(const Vector& x, const Vector& y) {
40 assert(x.size() == y.size());
41 double d = 0.0;
42 for (Int i = 0; i < (Int) x.size(); i++)
43 d += x[i]*y[i];
44 return d;
45 }
46
FindMaxAbs(const Vector & x)47 Int FindMaxAbs(const Vector& x) {
48 double xmax = 0.0;
49 Int imax = 0;
50 for (Int i = 0; i < (Int) x.size(); i++) {
51 if (std::abs(x[i]) > xmax) {
52 xmax = std::abs(x[i]);
53 imax = i;
54 }
55 }
56 return imax;
57 }
58
Permute(const std::vector<Int> & permuted_index,const Vector & rhs,Vector & lhs)59 void Permute(const std::vector<Int>& permuted_index, const Vector& rhs,
60 Vector& lhs) {
61 Int m = permuted_index.size();
62 for (Int i = 0; i < m; i++)
63 lhs[permuted_index[i]] = rhs[i];
64 }
65
PermuteBack(const std::vector<Int> & permuted_index,const Vector & rhs,Vector & lhs)66 void PermuteBack(const std::vector<Int>& permuted_index, const Vector& rhs,
67 Vector& lhs) {
68 Int m = permuted_index.size();
69 for (Int i = 0; i < m; i++)
70 lhs[i] = rhs[permuted_index[i]];
71 }
72
InversePerm(const std::vector<Int> & perm)73 std::vector<Int> InversePerm(const std::vector<Int>& perm) {
74 Int m = perm.size();
75 std::vector<Int> invperm(m);
76 // at() throws an exception if the index is out of range
77 for (Int i = 0; i < m; i++)
78 invperm.at(perm[i]) = i;
79 return invperm;
80 }
81
greater_or_equal(const std::pair<double,Int> & a,const std::pair<double,Int> & b)82 static bool greater_or_equal(const std::pair<double,Int>& a,
83 const std::pair<double,Int>& b) {
84 return a>b;
85 }
86
Sortperm(Int m,const double * values,bool reverse)87 std::vector<Int> Sortperm(Int m, const double* values, bool reverse) {
88 std::vector<Int> perm(m);
89 if (!values) {
90 for (Int i = 0; i < m; i++)
91 perm[i] = i;
92 return perm;
93 }
94 std::vector<std::pair<double,Int>> value_index(m);
95 for (Int i = 0; i < m; i++)
96 value_index[i] = std::make_pair(values[i], i);
97 if (reverse)
98 std::sort(value_index.begin(), value_index.end(), greater_or_equal);
99 else
100 std::sort(value_index.begin(), value_index.end());
101 for (Int i = 0; i < m; i++)
102 perm[i] = value_index[i].second;
103 return perm;
104 }
105
106 } // namespace ipx
107