1 /* linbox/field/modular-float.h
2  * Copyright (C) 2003 Pascal Giorgi
3  *               2007 Clement Pernet
4  * Written by Clement Pernet <cpernet@uwaterloo.ca>
5  *
6  * ------------------------------------
7  *
8  *
9  * ========LICENCE========
10  * This file is part of the library LinBox.
11  *
12  * LinBox is free software: you can redistribute it and/or modify
13  * it under the terms of the  GNU Lesser General Public
14  * License as published by the Free Software Foundation; either
15  * version 2.1 of the License, or (at your option) any later version.
16  *
17  * This library is distributed in the hope that it will be useful,
18  * but WITHOUT ANY WARRANTY; without even the implied warranty of
19  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
20  * Lesser General Public License for more details.
21  *
22  * You should have received a copy of the GNU Lesser General Public
23  * License along with this library; if not, write to the Free Software
24  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA
25  * ========LICENCE========
26  *.
27  */
28 
29 /*! @file field/modular/modular-float.h
30  * @ingroup field
31  * @brief  representation of <code>Z/mZ</code> over \c float .
32  */
33 
34 #ifndef __LINBOX_modular_float_H
35 #define __LINBOX_modular_float_H
36 
37 #ifdef __INTEL_COMPILER
38 #define FmodF fmodf
39 #else
40 #define FmodF fmod
41 #endif
42 
43 #include <cmath>
44 #include <givaro/modular-floating.h>
45 
46 #include "linbox/linbox-config.h"
47 #include "linbox/integer.h"
48 #include "linbox/vector/vector-domain.h"
49 #include "linbox/ring/modular.h"
50 #include "linbox/field/field-traits.h"
51 #include "linbox/util/field-axpy.h"
52 #include "linbox/util/debug.h"
53 #include "linbox/field/field-traits.h"
54 
55 // Namespace in which all LinBox code resides
56 namespace LinBox
57 {
58 	class MultiModFloat;
59 
60 	template <>
61 	class FieldAXPY<Givaro::Modular<float> > {
62 	public:
63 
64 		typedef float Element;
65 		typedef float Abnormal;
66 		typedef Givaro::Modular<float> Field;
67 
FieldAXPY(const Field & F)68 		FieldAXPY (const Field &F) :
69 			_field (&F) , //_invmod(1./field().fcharacteristic()),
70 			_y(0.) , _bound( (float) ( (1_i32 << 23) - (uint32_t) (field().characteristic()*field().characteristic())))
71 		{}
72 
FieldAXPY(const FieldAXPY & faxpy)73 		FieldAXPY (const FieldAXPY &faxpy) :
74 			_field (faxpy._field),// _invmod(faxpy._invmod) ,
75 			_y(faxpy._y), _bound(faxpy._bound)
76 		{}
77 
mulacc(const Element & a,const Element & x)78 		inline Element& mulacc (const Element &a, const Element &x)
79 		{
80 			Element tmp= a*x;
81 			return accumulate(tmp);
82 		}
83 
accumulate(const Element & tmp)84 		inline Element& accumulate (const Element &tmp)
85 		{
86 			_y += tmp;
87 			if (_y > _bound)
88 				return _y = fmodf (_y, field().fcharacteristic());
89 			else
90 				return _y;
91 		}
92 
get(Element & y)93 		inline Element& get (Element &y)
94 		{
95 			_y = fmodf (_y, field().fcharacteristic());
96 			return y=_y ;
97 		}
98 
assign(const Element y)99 		inline FieldAXPY &assign (const Element y)
100 		{
101 			_y = y;
102 			return *this;
103 		}
104 
reset()105 		inline void reset()
106 		{
107 			_y = 0.;
108 		}
109 
field()110 		inline const Field & field() const { return *_field; }
111 
112 	private:
113 
114 		const Field *_field;
115 		//float _invmod;
116 		float _y;
117 		float _bound;
118 	};
119 
120 
121 	template <>
122 	class DotProductDomain<Givaro::Modular<float> > : public VectorDomainBase<Givaro::Modular<float> > {
123 	private:
124 		float _bound;
125 		size_t _nmax;
126 		//float _invmod;
127 
128 	public:
129 		typedef float Element;
DotProductDomain(const Givaro::Modular<float> & F)130 		DotProductDomain (const Givaro::Modular<float> &F) :
131 			VectorDomainBase<Givaro::Modular<float> > (F)
132 			, _bound( (float) ( (1<<23) - (uint32_t) (F.characteristic()*F.characteristic())))
133 			//, _invmod(1./field().fcharacteristic())
134 		{
135 			_nmax= (size_t)floor((double(1<<11)* double(1<<12))/ double(F.fcharacteristic() * F.fcharacteristic()));
136 		}
137 
138 		using VectorDomainBase<Givaro::Modular<float> >::field;
139 	protected:
140 		template <class Vector1, class Vector2>
dotSpecializedDD(Element & res,const Vector1 & v1,const Vector2 & v2)141 		inline Element &dotSpecializedDD (Element &res, const Vector1 &v1, const Vector2 &v2) const
142 		{
143 
144 			float y = 0.;
145 			if (v1.size() < _nmax) {
146 				for (size_t i = 0; i< v1.size();++i)
147 					y += v1[i] * v2[i] ;
148 				y = fmodf(y, field().fcharacteristic());
149 			}
150 			else {
151 				float t = 0.;
152 				size_t i=0;
153 				for (;i< v1.size()- _nmax ;i=i+_nmax)
154 				{
155 					for (size_t j=i;j<i+_nmax;++j)
156 						y += v1[j] * v2[j];
157 					t+=fmodf(y, field().fcharacteristic());
158 					y=0.;
159 				}
160 				for (;i < v1.size();++i)
161 					y += v1[i] * v2[i];
162 				t+=fmodf(y, field().fcharacteristic());
163 				y = fmodf(t, field().fcharacteristic());
164 			}
165 			return res = y;
166 		}
167 
168 		template <class Vector1, class Vector2>
dotSpecializedDSP(Element & res,const Vector1 & v1,const Vector2 & v2)169 		inline Element &dotSpecializedDSP (Element &res, const Vector1 &v1, const Vector2 &v2) const
170 		{
171 
172 			float y = 0.;
173 
174 
175 			if (v1.first.size() < _nmax)
176 			{
177 				for (size_t i=0;i<v1.first.size();++i)
178 					y+= v1.second[i] * v2[v1.first[i]];
179 				y = fmodf(y, field().fcharacteristic());
180 			}
181 			else
182 			{
183 				float t =0.;
184 				size_t i=0;
185 				for (;i< v1.first.size()- _nmax ;i=i+_nmax)
186 				{
187 					for (size_t j=i;j<i+_nmax;++j)
188 						y += v1.second[j] * v2[v1.first[j]];
189 					t+=fmodf(y, field().fcharacteristic());
190 					y=0.;
191 				}
192 				for (;i < v1.first.size();++i)
193 					y += v1.second[i] * v2[v1.first[i]];
194 				t+= fmodf(y, field().fcharacteristic());
195 				y = fmodf(t, field().fcharacteristic());
196 			}
197 			return res = y;
198 		}
199 	};
200 }
201 
202 #undef FmodF
203 
204 #endif //__LINBOX_modular_float_H
205 
206 // Local Variables:
207 // mode: C++
208 // tab-width: 4
209 // indent-tabs-mode: nil
210 // c-basic-offset: 4
211 // End:
212 // vim:sts=4:sw=4:ts=4:et:sr:cino=>s,f0,{0,g0,(0,\:0,t0,+0,=s
213