File: permutation.cweb

package info (click to toggle)
dynare 4.5.7-1
  • links: PTS, VCS
  • area: main
  • in suites: buster
  • size: 49,408 kB
  • sloc: cpp: 84,998; ansic: 29,058; pascal: 13,843; sh: 4,833; objc: 4,236; yacc: 3,622; makefile: 2,278; lex: 1,541; python: 236; lisp: 69; xml: 8
file content (188 lines) | stat: -rw-r--r-- 4,553 bytes parent folder | download | duplicates (5)
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
@q $Id: permutation.cweb 332 2005-07-15 13:41:48Z kamenik $ @>
@q Copyright 2004, Ondra Kamenik @>

@ Start of {\tt permutation.cweb} file.
@c

#include "permutation.h"
#include "tl_exception.h"

@<|Permutation::apply| code@>;
@<|Permutation::inverse| code@>;
@<|Permutation::tailIdentity| code@>;
@<|Permutation::computeSortingMap| code@>;
@<|PermutationSet| constructor code 1@>;
@<|PermutationSet| constructor code 2@>;
@<|PermutationSet| destructor code@>;
@<|PermutationSet::getPreserving| code@>;
@<|PermutationBundle| constructor code@>;
@<|PermutationBundle| destructor code@>;
@<|PermutationBundle::get| code@>;
@<|PermutationBundle::generateUpTo| code@>;


@ This is easy, we simply apply the map in the fashion $s\circ m$..
@<|Permutation::apply| code@>=
void Permutation::apply(const IntSequence& src, IntSequence& tar) const
{
	TL_RAISE_IF(src.size() != permap.size() || tar.size() != permap.size(),
				"Wrong sizes of input or output in Permutation::apply");
	for (int i = 0; i < permap.size(); i++)
		tar[i] = src[permap[i]];
}


void Permutation::apply(IntSequence& tar) const
{
	IntSequence tmp(tar);
	apply(tmp, tar);
}

@ 
@<|Permutation::inverse| code@>=
void Permutation::inverse()
{
	IntSequence former(permap);
	for (int i = 0; i < size(); i++)
		permap[former[i]] = i;
}


@ Here we find a number of trailing indices which are identical with
the permutation.

@<|Permutation::tailIdentity| code@>=
int Permutation::tailIdentity() const
{
	int i = permap.size();
	while (i > 0 && permap[i-1] == i-1)
		i--;
	return permap.size() - i;
}

@ This calculates a map which corresponds to sorting in the following
sense: $(\hbox{sorted }s)\circ m = s$, where $s$ is a given sequence.

We go through |s| and find an the same item in sorted |s|. We
construct the |permap| from the found pair of indices. We have to be
careful, to not assign to two positions in |s| the same position in
sorted |s|, so we maintain a bitmap |flag|, in which we remember
indices from the sorted |s| already assigned.

@<|Permutation::computeSortingMap| code@>=
void Permutation::computeSortingMap(const IntSequence& s)
{
	IntSequence srt(s);
	srt.sort();
	IntSequence flags(s.size(),0);

	for (int i = 0; i < s.size(); i++) {
		int j = 0;
		while (j < s.size() && (flags[j] || srt[j] != s[i]))
			j++;
		TL_RAISE_IF(j == s.size(),
					"Internal algorithm error in Permutation::computeSortingMap");
		flags[j] = 1;
		permap[i] = j;
	}
}

@ 
@<|PermutationSet| constructor code 1@>=
PermutationSet::PermutationSet()
	: order(1), size(1), pers(new const Permutation*[size])
{
	pers[0] = new Permutation(1);
}

@ 
@<|PermutationSet| constructor code 2@>=
PermutationSet::PermutationSet(const PermutationSet& sp, int n)
	: order(n), size(n*sp.size),
	  pers(new const Permutation*[size])
{
	for (int i = 0; i < size; i++)
		pers[i] = NULL;

	TL_RAISE_IF(n != sp.order+1,
				"Wrong new order in PermutationSet constructor");

	int k = 0;
	for (int i = 0; i < sp.size; i++) {
		for (int j = 0;	 j < order; j++,k++) {
			pers[k] = new Permutation(*(sp.pers[i]), j);
		}
	}
}

@ 
@<|PermutationSet| destructor code@>=
PermutationSet::~PermutationSet()
{
	for (int i = 0; i < size; i++)
		if (pers[i])
			delete pers[i];
	delete [] pers;
}

@ 
@<|PermutationSet::getPreserving| code@>=
vector<const Permutation*> PermutationSet::getPreserving(const IntSequence& s) const
{
	TL_RAISE_IF(s.size() != order,
				"Wrong sequence length in PermutationSet::getPreserving");

	vector<const Permutation*> res;
	IntSequence tmp(s.size());
	for (int i = 0; i < size; i++) {
		pers[i]->apply(s, tmp);
		if (s == tmp) {
			res.push_back(pers[i]);
		}
	}

	return res;
}

@ 
@<|PermutationBundle| constructor code@>=
PermutationBundle::PermutationBundle(int nmax)
{
	nmax = max(nmax, 1);
	generateUpTo(nmax);
}

@ 
@<|PermutationBundle| destructor code@>=
PermutationBundle::~PermutationBundle()
{
	for (unsigned int i = 0; i < bundle.size(); i++)
		delete bundle[i];
}

@ 
@<|PermutationBundle::get| code@>=
const PermutationSet& PermutationBundle::get(int n) const
{
	if (n > (int)(bundle.size()) || n < 1) {
		TL_RAISE("Permutation set not found in PermutationSet::get");
		return *(bundle[0]);
	} else {
		return *(bundle[n-1]);
	}
}

@ 
@<|PermutationBundle::generateUpTo| code@>=
void PermutationBundle::generateUpTo(int nmax)
{
	if (bundle.size() == 0)
		bundle.push_back(new PermutationSet());

	int curmax = bundle.size();
	for (int n = curmax+1; n <= nmax; n++) {
		bundle.push_back(new PermutationSet(*(bundle.back()), n));
	}
}

@ End of {\tt permutation.cweb} file.