1// ========================================================================== 2// Copyright(c)'1994-2015 by The Givaro group 3// This file is part of Givaro. 4// Givaro is governed by the CeCILL-B license under French law 5// and abiding by the rules of distribution of free software. 6// see the COPYRIGHT file for more details. 7// Authors: Clement Pernet <clement.pernet@gmail.com> 8// A. Breust (taken from FFLAS-FFPACK) 9// ========================================================================== 10 11#ifndef __GIVARO_modular_balanced_float_INL 12#define __GIVARO_modular_balanced_float_INL 13 14#include <cmath> // fmod 15 16#define NORMALISE(x) \ 17{ \ 18 if (x < _mhalfp) x += _p; \ 19 else if (x > _halfp) x -= _p; \ 20} 21 22#define NORMALISE_HI(x) \ 23{ \ 24 if (x > _halfp) x -= _p; \ 25} 26 27namespace Givaro 28{ 29 30 //----- Classic arithmetic 31 32 inline ModularBalanced<float>::Element& 33 ModularBalanced<float>::mul(Element& r, const Element& a, const Element& b) const 34 { 35 return reduce(r = a * b); 36 } 37 38 inline ModularBalanced<float>::Element& 39 ModularBalanced<float>::div(Element& r, const Element& a, const Element& b) const 40 { 41 Element tmp; 42 return mul (r, a, inv(tmp, b)); 43 } 44 45 inline ModularBalanced<float>::Element& 46 ModularBalanced<float>::add(Element& r, const Element& a, const Element& b) const 47 { 48 r = a + b; 49 NORMALISE(r); 50 return r; 51 } 52 53 inline ModularBalanced<float>::Element& 54 ModularBalanced<float>::sub(Element& r, const Element& a, const Element& b) const 55 { 56 r = a - b; 57 NORMALISE(r); 58 return r; 59 } 60 61 inline ModularBalanced<float>::Element& 62 ModularBalanced<float>::neg(Element& r, const Element& a) const 63 { 64 return r = -a; 65 } 66 67 inline ModularBalanced<float>::Element& 68 ModularBalanced<float>::inv(Element& r, const Element& a) const 69 { 70 r = invext(a, _p); 71 NORMALISE(r); 72 return r; 73 } 74 75 inline bool ModularBalanced<float>::isUnit(const Element& a) const 76 { 77 Element u,d; 78 extended_euclid(u,d,a,_p); 79 return isOne(d) || isMOne(d); 80 } 81 82 inline ModularBalanced<float>::Element& 83 ModularBalanced<float>::mulin(Element& r, const Element& a) const 84 { 85 return mul(r, r, a); 86 } 87 88 inline ModularBalanced<float>::Element& 89 ModularBalanced<float>::divin(Element& r, const Element& a) const 90 { 91 return div(r, r, a); 92 } 93 94 inline ModularBalanced<float>::Element& 95 ModularBalanced<float>::addin(Element& r, const Element& a) const 96 { 97 return add(r, r, a); 98 } 99 100 inline ModularBalanced<float>::Element& 101 ModularBalanced<float>::subin(Element& r, const Element& a) const 102 { 103 return sub(r, r, a); 104 } 105 106 inline ModularBalanced<float>::Element& 107 ModularBalanced<float>::negin(Element& r) const 108 { 109 return neg(r, r); 110 } 111 112 inline ModularBalanced<float>::Element& 113 ModularBalanced<float>::invin(Element& r) const 114 { 115 return inv(r, r); 116 } 117 118 //----- Special arithmetic 119 120 inline ModularBalanced<float>::Element& 121 ModularBalanced<float>::axpy(Element& r, const Element& a, const Element& x, const Element& y) const 122 { 123 return reduce(r = a * x + y); 124 } 125 126 inline ModularBalanced<float>::Element& 127 ModularBalanced<float>::axpyin(Element& r, const Element& a, const Element& x) const 128 { 129 return reduce(r += a * x); 130 } 131 132 inline ModularBalanced<float>::Element& 133 ModularBalanced<float>::axmy(Element& r, const Element& a, const Element& x, const Element& y) const 134 { 135 return reduce(r = a * x - y); 136 } 137 138 inline ModularBalanced<float>::Element& 139 ModularBalanced<float>::axmyin(Element& r, const Element& a, const Element& x) const 140 { 141 return reduce(r = a * x - r); 142 } 143 144 inline ModularBalanced<float>::Element& 145 ModularBalanced<float>::maxpy(Element& r, const Element& a, const Element& x, const Element& y) const 146 { 147 return reduce(r = y - a * x); 148 } 149 150 inline ModularBalanced<float>::Element& 151 ModularBalanced<float>:: maxpyin(Element& r, const Element& a, const Element& x) const 152 { 153 return reduce(r -= a * x); 154 } 155 156 //----- Init 157 158 inline ModularBalanced<float>::Element& 159 ModularBalanced<float>::init(Element& x) const 160 { 161 return x = 0; 162 } 163 164 inline ModularBalanced<float>::Element& 165 ModularBalanced<float>::init(Element& x, const float y) const 166 { 167 x = static_cast<Element>(fmodf(y, _p)); 168 NORMALISE(x); 169 return x; 170 } 171 172 inline ModularBalanced<float>::Element& 173 ModularBalanced<float>::init(Element& x, const double y) const 174 { 175 x = static_cast<Element>(fmod(y, double(_p))); 176 NORMALISE(x); 177 return x; 178 } 179 180 inline ModularBalanced<float>::Element& 181 ModularBalanced<float>::init(Element& x, const int32_t y) const 182 { 183 x = static_cast<Element>(y % static_cast<int32_t>(_up)); 184 NORMALISE(x); 185 return x; 186 } 187 188 inline ModularBalanced<float>::Element& 189 ModularBalanced<float>::init(Element& x, const uint32_t y) const 190 { 191 x = static_cast<Element>(y % static_cast<uint32_t>(_up)); 192 NORMALISE_HI(x); 193 return x; 194 } 195 196 inline ModularBalanced<float>::Element& 197 ModularBalanced<float>::init(Element& x, const int64_t y) const 198 { 199 x = static_cast<Element>(y % static_cast<int64_t>(_up)); 200 NORMALISE(x); 201 return x; 202 } 203 204 inline ModularBalanced<float>::Element& 205 ModularBalanced<float>::init(Element& x, const uint64_t y) const 206 { 207 x = static_cast<Element>(y % static_cast<uint64_t>(_up)); 208 NORMALISE_HI(x); 209 return x; 210 } 211 212 inline ModularBalanced<float>::Element& 213 ModularBalanced<float>::init(Element& x, const Integer& y) const 214 { 215 x = static_cast<Element>(y % _p); 216 NORMALISE_HI(x); 217 return x; 218 } 219 220 inline ModularBalanced<float>::Element& 221 ModularBalanced<float>::assign(Element& x, const Element& y) const 222 { 223 return x = y; 224 } 225 226 //----- Reduce 227 228 inline ModularBalanced<float>::Element& 229 ModularBalanced<float>::reduce(Element& x, const Element& y) const 230 { 231 x = fmodf(y, _p); 232 NORMALISE(x); 233 return x; 234 } 235 236 inline ModularBalanced<float>::Element& 237 ModularBalanced<float>::reduce(Element& x) const 238 { 239 x = fmodf(x, _p); 240 NORMALISE(x); 241 return x; 242 } 243 244 //----- IO 245 246 inline std::ostream& 247 ModularBalanced<float>::write(std::ostream& os) const 248 { 249 return os << "ModularBalanced<float> mod " << _p; 250 } 251 252 inline std::ostream& 253 ModularBalanced<float>::write(std::ostream& os, const Element& x) const 254 { 255 return os << x; 256 } 257 258 inline std::istream& 259 ModularBalanced<float>::read(std::istream& is, Element& x) const 260 { 261 Element tmp; 262 is >> tmp; 263 init(x, tmp); 264 return is; 265 } 266 267} 268 269#undef NORMALISE 270#undef NORMALISE_HI 271 272#endif 273/* -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */ 274// vim:sts=4:sw=4:ts=4:et:sr:cino=>s,f0,{0,g0,(0,\:0,t0,+0,=s 275