This file is indexed.

/usr/include/linbox/algorithms/cra-full-multip.h is in liblinbox-dev 1.1.6~rc0-4.1.

This file is owned by root:root, with mode 0o644.

The actual contents of the file can be viewed below.

  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
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
/* -*- mode: C++; tab-width: 8; indent-tabs-mode: t; c-basic-offset: 8 -*- */
// ======================================================================= //
// Time-stamp: <23 Apr 08 19:28:37 Jean-Guillaume.Dumas@imag.fr> 
// ======================================================================= //
#ifndef __LINBOX_CRA_FULL_MULTIP_H
#define __LINBOX_CRA_FULL_MULTIP_H

#include "linbox/util/timer.h"
#include <stdlib.h>
#include "linbox/integer.h"
#include "linbox/solutions/methods.h"
#include <vector>
#include <utility>

#include "linbox/algorithms/lazy-product.h"

namespace LinBox {
    
template<class Domain_Type>
struct FullMultipCRA {
	typedef Domain_Type			Domain;
	typedef typename Domain::Element DomainElement;
	typedef FullMultipCRA<Domain> 		Self_t;
	
protected:
	std::vector< double >           	RadixSizes_;
	std::vector< LazyProduct >      	RadixPrimeProd_;
	std::vector< std::vector<Integer> >    	RadixResidues_;
	std::vector< bool >             	RadixOccupancy_;
	double                			LOGARITHMIC_UPPER_BOUND;
	
	
public:
	FullMultipCRA(const double b=0.0) : LOGARITHMIC_UPPER_BOUND(b) {}
        
	template< template<class, class> class Vect, template <class> class Alloc>
	void initialize (const Domain& D, const Vect<DomainElement, Alloc<DomainElement> >& e) {
		RadixSizes_.resize(1);
		RadixPrimeProd_.resize(1);
		RadixResidues_.resize(1);
		RadixOccupancy_.resize(1); RadixOccupancy_.front() = false;
		progress(D, e);
	}
        
	template< template<class, class> class Vect, template <class> class Alloc>
	void progress (const Domain& D, const Vect<DomainElement, Alloc<DomainElement> >& e) {
		// Radix shelves
		std::vector< double >::iterator  _dsz_it = RadixSizes_.begin();
		std::vector< LazyProduct >::iterator _mod_it = RadixPrimeProd_.begin();
		std::vector< std::vector<Integer> >::iterator _tab_it = RadixResidues_.begin();
		std::vector< bool >::iterator    _occ_it = RadixOccupancy_.begin();
		std::vector<Integer> ri(e.size()); LazyProduct mi; double di;
		if (*_occ_it) {
			// If lower shelf is occupied
			// Combine it with the new residue
			// The for loop will try to put the resulting combination on the upper shelf
			typename Vect<DomainElement, Alloc<DomainElement> >::const_iterator  e_it = e.begin();
			std::vector<Integer>::iterator       ri_it = ri.begin();
			std::vector<Integer>::const_iterator t0_it = _tab_it->begin();
			for( ; e_it != e.end(); ++e_it, ++ri_it, ++ t0_it)
				fieldreconstruct(*ri_it, D, *e_it, *t0_it, (*_mod_it).operator()() );
			Integer tmp; D.characteristic(tmp);
			di = *_dsz_it + log(double(tmp)); 
			mi.mulin(tmp);
			mi.mulin(*_mod_it);
			*_occ_it = false;
		} else {
			// Lower shelf is free
			// Put the new residue here and exit
			Integer tmp; D.characteristic(tmp);
			_mod_it->initialize(tmp);
			*_dsz_it = log(double(tmp));
			typename Vect<DomainElement, Alloc<DomainElement> >::const_iterator e_it = e.begin();
			_tab_it->resize(e.size());
			std::vector<Integer>::iterator t0_it= _tab_it->begin();
			for( ; e_it != e.end(); ++e_it, ++ t0_it)
				D.convert(*t0_it, *e_it);
			*_occ_it = true;
			return;
		}
		
		// We have a combination to put in the upper shelf
		for(++_dsz_it, ++_mod_it, ++_tab_it, ++_occ_it ; _occ_it != RadixOccupancy_.end() ; ++_dsz_it, ++_mod_it, ++_tab_it, ++_occ_it) {
			if (*_occ_it) {
				// This shelf is occupied
				// Combine it with the new combination
				// The loop will try to put it on the upper shelf
				std::vector<Integer>::iterator      ri_it = ri.begin();
				std::vector<Integer>::const_iterator t_it= _tab_it->begin();
                                
                                Integer invprod; precomputeInvProd(invprod, mi(), _mod_it->operator()());
				for( ; ri_it != ri.end(); ++ri_it, ++ t_it) 
                                    smallbigreconstruct(*ri_it, *t_it, invprod); 

                                    // Product (lazy) computation
				mi.mulin(*_mod_it);

                                    // Moding out
                                for(ri_it = ri.begin() ; ri_it != ri.end(); ++ri_it) {
                                    *ri_it %= mi();
                                }

				di += *_dsz_it;
				*_occ_it = false;
			} else {
				// This shelf is free
				// Put the new combination here and exit
				*_dsz_it = di;
				*_mod_it = mi;
				*_tab_it = ri;
				*_occ_it = true;
				return;
			}
		}
		// All the shelfves were occupied
		// We create a new top shelf 
		// And put the new combination there
		RadixSizes_.push_back( di );
		RadixResidues_.push_back( ri );
		RadixPrimeProd_.push_back( mi );
		RadixOccupancy_.push_back ( true );
	}
	
        
	
	template<template<class, class> class Vect, template <class> class Alloc>
	Vect<Integer, Alloc<Integer> >& result (Vect<Integer, Alloc<Integer> > &d){
		d.resize( (RadixResidues_.front()).size() );
		std::vector< LazyProduct >::iterator          _mod_it = RadixPrimeProd_.begin();
		std::vector< std::vector< Integer > >::iterator _tab_it = RadixResidues_.begin();
		std::vector< bool >::iterator                _occ_it = RadixOccupancy_.begin();
		LazyProduct Product;
		// We have to find to lowest occupied shelf
		for( ; _occ_it != RadixOccupancy_.end() ; ++_mod_it, ++_tab_it, ++_occ_it) {
			if (*_occ_it) {
				// Found the lowest occupied shelf
				Product = *_mod_it;
				std::vector<Integer>::iterator t0_it = d.begin();
				std::vector<Integer>::iterator t_it = _tab_it->begin();
				if (++_occ_it == RadixOccupancy_.end()) {
					// It is the only shelf of the radix 
					// We normalize the result and output it
					for( ; t0_it != d.end(); ++t0_it, ++t_it)
						normalize(*t0_it = *t_it, *t_it, _mod_it->operator()());
					RadixPrimeProd_.resize(1);
					return d;
				} else {
					// There are other shelves
					// The result is initialized with this shelf
					// The for loop will combine the other shelves m with the actual one
					for( ; t0_it != d.end(); ++t0_it, ++t_it)
						*t0_it  = *t_it;
					++_mod_it; ++_tab_it; 
					break;
				}
			}
		}
		for( ; _occ_it != RadixOccupancy_.end() ; ++_mod_it, ++_tab_it, ++_occ_it) {
			if (*_occ_it) {
				// This shelf is occupied
				// We need to combine it with the actual value of the result
				std::vector<Integer>::iterator t0_it = d.begin();
				std::vector<Integer>::const_iterator t_it = _tab_it->begin();
                                Integer invprod; 
                                precomputeInvProd(invprod, Product(), _mod_it->operator()() );

				for( ; t0_it != d.end(); ++t0_it, ++t_it)
                                    smallbigreconstruct(*t0_it, *t_it, invprod);
                                
                                    // Overall product computation
				Product.mulin(*_mod_it);

                                    // Moding out and normalization
                                for(t0_it = d.begin();t0_it != d.end(); ++t0_it) {
                                    *t0_it %= Product();
                                    Integer tmp(*t0_it);
                                    normalize(*t0_it, tmp, Product());
                                }
                                
			}
		}

		// We put it also the final prime product in the first shelf of products
		// JGD : should we also put the result 
		//       in the first shelf of residues and resize it to 1
		//       and set to true the first occupancy and resize it to 1
		//       in case result is not the last call (more progress to go) ?
		RadixPrimeProd_.resize(1);
		RadixPrimeProd_.front() = Product;
		return d;
	}
	
	bool terminated() {
		double logm(0.0);
		std::vector< double >::iterator _dsz_it = RadixSizes_.begin();
		std::vector< bool >::iterator _occ_it = RadixOccupancy_.begin();
		for( ; _occ_it != RadixOccupancy_.end() ; ++_dsz_it, ++_occ_it) {
			if (*_occ_it) {
				logm = *_dsz_it;
				++_dsz_it; ++_occ_it;
				break;
			}
		}
		for( ; _occ_it != RadixOccupancy_.end() ; ++_dsz_it, ++_occ_it) {
			if (*_occ_it) {
				logm += *_dsz_it;
			}
		}
		return logm > LOGARITHMIC_UPPER_BOUND;
	}
	
	
	
	bool noncoprime(const Integer& i) const {
		std::vector< LazyProduct >::const_iterator _mod_it = RadixPrimeProd_.begin();
		std::vector< bool >::const_iterator    _occ_it = RadixOccupancy_.begin();
		for( ; _occ_it != RadixOccupancy_.end() ; ++_mod_it, ++_occ_it)
			if ((*_occ_it) && (_mod_it->noncoprime(i))) return true;
		return false;
	}
	
	

	
protected:
	
    Integer& precomputeInvProd(Integer& res, const Integer& m1, const Integer& m0) {
        inv(res, m0, m1);
        return res *= m0; // res <-- (m0^{-1} mod m1) m0
    }

    Integer& smallbigreconstruct(Integer& u1, const Integer& u0, const Integer& invprod) {
        u1 -= u0;	  // u1 <-- (u1-u0)
        u1 *= invprod;    // u1 <-- (u1-u0)( m0^{-1} mod m1 ) m0 
        return u1 += u0;  // u1 <-- u0 + (u1-u0)( m0^{-1} mod m1 ) m0
    }
    
        
	Integer& normalize(Integer& u1, Integer& tmp, const Integer& m1) {
		if (u1 < 0)
			tmp += m1;
		else
			tmp -= m1;
		return ((absCompare(u1,tmp) > 0)? u1 = tmp : u1 );
	}
	
	
	Integer& fieldreconstruct(Integer& res, const Domain& D1, const DomainElement& u1, const Integer& r0, const Integer& P0) {
		DomainElement u0, m0; D1.init(u0, r0);
		if (D1.areEqual(u1, u0))
			return res=r0;
		else
			return fieldreconstruct(res, D1, u1, u0, D1.init(m0, P0), r0, P0);
	}
	
	Integer& fieldreconstruct(Integer& res, const Domain& D1, const DomainElement& u1, DomainElement& u0, DomainElement& m0, const Integer& r0, const Integer& P0) {
		// u0 and m0 are modified
		D1.negin(u0);   	// u0 <-- -u0
		D1.addin(u0,u1);   	// u0 <-- u1-u0
		D1.invin(m0);   	// m0 <-- m0^{-1} mod m1
		D1.mulin(u0, m0);   // u0 <-- (u1-u0)( m0^{-1} mod m1 )
		D1.convert(res, u0);// res <-- (u1-u0)( m0^{-1} mod m1 )         and res <  m1
		res *= P0;      	// res <-- (u1-u0)( m0^{-1} mod m1 ) m0      and res <= (m0m1-m0)
		return res += r0;	// res <-- u0 + (u1-u0)( m0^{-1} mod m1 ) m0 and res <  m0m1
	}
	
	
};
	
}


#endif