1 /* 2 * Copyright © 2004 Ondra Kamenik 3 * Copyright © 2019 Dynare Team 4 * 5 * This file is part of Dynare. 6 * 7 * Dynare is free software: you can redistribute it and/or modify 8 * it under the terms of the GNU General Public License as published by 9 * the Free Software Foundation, either version 3 of the License, or 10 * (at your option) any later version. 11 * 12 * Dynare is distributed in the hope that it will be useful, 13 * but WITHOUT ANY WARRANTY; without even the implied warranty of 14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 15 * GNU General Public License for more details. 16 * 17 * You should have received a copy of the GNU General Public License 18 * along with Dynare. If not, see <http://www.gnu.org/licenses/>. 19 */ 20 21 // Equivalences. 22 23 /* Here we define an equivalence of a set of integers {0, 1, …, k-1}. 24 The purpose is clear, in the tensor library we often iterate 25 through all equivalences and sum matrices. We need an abstraction for 26 an equivalence class, equivalence and a set of all equivalences. 27 28 The equivalence class (which is basically a set of integers) is here 29 implemented as ordered integer sequence. The ordered sequence is not 30 implemented via IntSequence, but via vector<int> since we need 31 insertions. The equivalence is implemented as an ordered list of 32 equivalence classes, and equivalence set is a list of equivalences. 33 34 The ordering of the equivalence classes within an equivalence is very 35 important. For instance, if we iterate through equivalences for k=5 36 and pickup some equivalence class, say { {0,4}, {1,2}, {3} }, we 37 then evaluate something like: 38 39 [B_y²u³]_α₁α₂β₁β₂β₃ = … + [f_z³]_γ₁γ₂γ₃ [z_yu]^γ₁_α₁β₃ [z_yu]^γ₂_α₂β₁ [zᵤ]^γ₃_β₂ + … 40 41 If the tensors are unfolded, we can evaluate this expression as 42 f_z³·(z_yu ⊗ z_yu ⊗ zᵤ)·P, where P is a suitable permutation of columns of 43 the expressions, which permutes them so that the index 44 (α₁,β₃,α₂,β₁,β₂) would go to (α₁,α₂,β₁,β₂,β₃). 45 46 The permutation P can be very ineffective (copying great amount of 47 small chunks of data) if the equivalence class ordering is chosen 48 badly. However, we do not provide any heuristic minimizing a total 49 time spent in all permutations. We choose an ordering which orders the 50 classes according to their averages, and according to the smallest 51 equivalence class element if the averages are the same. */ 52 53 #ifndef EQUIVALENCE_H 54 #define EQUIVALENCE_H 55 56 #include "int_sequence.hh" 57 58 #include <vector> 59 #include <list> 60 #include <string> 61 62 /* Here is the abstraction for an equivalence class. We implement it as 63 vector<int>. We have a constructor for empty class, copy 64 constructor. What is important here is the ordering operator 65 operator<() and methods for addition of an integer, and addition of 66 another sequence. Also we provide method has() which returns true if a 67 given integer is contained. */ 68 69 class OrdSequence 70 { 71 std::vector<int> data; 72 public: OrdSequence()73 OrdSequence() : data() 74 { 75 } 76 bool operator==(const OrdSequence &s) const; 77 int operator[](int i) const; 78 bool operator<(const OrdSequence &s) const; 79 const std::vector<int> & getData() const80 getData() const 81 { 82 return data; 83 } 84 int length() const85 length() const 86 { 87 return data.size(); 88 } 89 void add(int i); 90 void add(const OrdSequence &s); 91 bool has(int i) const; 92 void print(const std::string &prefix) const; 93 private: 94 double average() const; 95 }; 96 97 /* Here is the abstraction for the equivalence. It is a list of 98 equivalence classes. Also we remember n, which is a size of 99 underlying set {0, 1, …, n-1}. 100 101 Method trace() “prints” the equivalence into the integer sequence. */ 102 103 class Permutation; 104 class Equivalence 105 { 106 private: 107 int n; 108 std::list<OrdSequence> classes; 109 public: 110 using const_seqit = std::list<OrdSequence>::const_iterator; 111 using seqit = std::list<OrdSequence>::iterator; 112 113 // Constructs { {0}, {1}, …, {n-1} } 114 explicit Equivalence(int num); 115 // Constructs { {0,1,…,n-1 } } 116 Equivalence(int num, const std::string &dummy); 117 // Copy constructor plus gluing i1 and i2 in one class 118 Equivalence(const Equivalence &e, int i1, int i2); 119 120 bool operator==(const Equivalence &e) const; 121 bool operator !=(const Equivalence & e) const122 operator!=(const Equivalence &e) const 123 { 124 return !operator==(e); 125 } 126 int getN() const127 getN() const 128 { 129 return n; 130 } 131 int numClasses() const132 numClasses() const 133 { 134 return classes.size(); 135 } 136 void trace(IntSequence &out, int n) const; 137 void trace(IntSequence & out) const138 trace(IntSequence &out) const 139 { 140 trace(out, numClasses()); 141 } 142 void trace(IntSequence &out, const Permutation &per) const; 143 void print(const std::string &prefix) const; 144 seqit begin()145 begin() 146 { 147 return classes.begin(); 148 } 149 const_seqit begin() const150 begin() const 151 { 152 return classes.begin(); 153 } 154 seqit end()155 end() 156 { 157 return classes.end(); 158 } 159 const_seqit end() const160 end() const 161 { 162 return classes.end(); 163 } 164 const_seqit find(int i) const; 165 seqit find(int i); 166 protected: 167 /* Here we have find methods. We can find an equivalence class having a 168 given number or we can find an equivalence class of a given index within 169 the ordering. 170 171 We have also an insert() method which inserts a given class 172 according to the class ordering. */ 173 const_seqit findHaving(int i) const; 174 seqit findHaving(int i); 175 void insert(const OrdSequence &s); 176 177 }; 178 179 /* The EquivalenceSet is a list of equivalences. The unique 180 constructor constructs a set of all equivalences over an n-elements 181 set. The equivalences are sorted in the list so that equivalences with 182 fewer number of classes are in the end. 183 184 The two methods has() and addParents() are useful in the constructor. */ 185 186 class EquivalenceSet 187 { 188 int n; 189 std::list<Equivalence> equis; 190 public: 191 explicit EquivalenceSet(int num); 192 void print(const std::string &prefix) const; 193 auto begin() const194 begin() const 195 { 196 return equis.begin(); 197 } 198 auto end() const199 end() const 200 { 201 return equis.end(); 202 } 203 private: 204 bool has(const Equivalence &e) const; 205 void addParents(const Equivalence &e, std::list<Equivalence> &added); 206 }; 207 208 /* The equivalence bundle class only encapsulates EquivalenceSet·s 209 from 1 up to a given number. It is able to retrieve the equivalence set 210 over n-element set for a given n, and also it can generate some more 211 sets on request. 212 213 It is fully responsible for storage needed for EquivalenceSet·s. */ 214 215 class EquivalenceBundle 216 { 217 std::vector<EquivalenceSet> bundle; 218 public: 219 explicit EquivalenceBundle(int nmax); 220 const EquivalenceSet &get(int n) const; 221 void generateUpTo(int nmax); 222 }; 223 224 #endif 225