PageRenderTime 52ms CodeModel.GetById 24ms RepoModel.GetById 0ms app.codeStats 0ms

/givaro-3.7.0/src/library/vector/givvectdenseops.inl

#
C++ Header | 441 lines | 335 code | 58 blank | 48 comment | 90 complexity | 329077f75ada8bdae6830854894be70e MD5 | raw file
Possible License(s): LGPL-3.0
  1. // ==========================================================================
  2. // $Source
  3. // Copyright(c)'1994-2009 by The Givaro group
  4. // This file is part of Givaro.
  5. // Givaro is governed by the CeCILL-B license under French law
  6. // and abiding by the rules of distribution of free software.
  7. // see the COPYRIGHT file for more details.
  8. // Authors: T. Gautier
  9. // $Id
  10. // ==========================================================================
  11. namespace Givaro {
  12. #pragma message "#warning this file will probably not compile"
  13. // -- map of a unary operator, with operator()( Type_t& res )
  14. // res and u could be aliases if OP permits it
  15. template<class Domain>
  16. template<class UNOP>
  17. inline void VectorDom<Domain,Dense>::
  18. map ( Rep& res, UNOP& op ) const
  19. {
  20. size_t sz = dim(res);
  21. for (size_t i=0; i<sz; ++i) op(res[i]);
  22. }
  23. template<class Domain>
  24. template<class BINOP>
  25. inline void VectorDom<Domain,Dense>::
  26. map( Rep& res, const BINOP& OP, const Rep& u, const Rep& v ) const
  27. {
  28. GIVARO_ASSERT((dim(res) == dim(u)) && (dim(u) == dim(v)), "Bad size");
  29. size_t sz = dim(res);
  30. for (size_t i=0; i<sz; ++i)
  31. OP(res[i], u[i], v[i]);
  32. }
  33. template<class Domain>
  34. template<class UNOP>
  35. inline void VectorDom<Domain,Dense>::
  36. map( Rep& res, UNOP& OP, const Rep& u ) const
  37. {
  38. GIVARO_ASSERT((dim(res) == dim(u)), "Bad size");
  39. size_t sz = dim(u);
  40. for (size_t i=0; i<sz; ++i)
  41. OP(res[i], u[i]);
  42. }
  43. // --- mul
  44. template<class Domain>
  45. inline void VectorDom<Domain,Dense>::mulin( Rep& res, const Type_t& val ) const
  46. {
  47. GIVARO_ASSERT((dim(res) == dim(u)), "Bad size");
  48. Curried2<MulOp<Domain> > opcode(_domain, (Type_t&)val);
  49. map( res, opcode);
  50. }
  51. template<class Domain>
  52. inline void VectorDom<Domain,Dense>::mul
  53. ( Rep& res, const Rep& u, const Type_t& val ) const
  54. {
  55. GIVARO_ASSERT((dim(res) == dim(u)), "Bad size");
  56. Curried2<MulOp<Domain> > opcode(_domain, (Type_t&)val);
  57. map( res, opcode, u);
  58. }
  59. template<class Domain>
  60. inline void VectorDom<Domain,Dense>::mul
  61. ( Rep& res, const Type_t& val, const Rep& v ) const
  62. {
  63. GIVARO_ASSERT((dim(res) == dim(v)), "Bad size");
  64. Curried1<MulOp<Domain> > opcode(_domain, val);
  65. map( res, opcode, v);
  66. }
  67. // --- add
  68. template<class Domain>
  69. inline void VectorDom<Domain,Dense>::addin( Rep& res, const Rep& u ) const
  70. {
  71. GIVARO_ASSERT((dim(res) == dim(u)), "Bad size");
  72. AddOp<Domain> opcode( _domain );
  73. map( res, opcode, u );
  74. }
  75. template<class Domain>
  76. inline void VectorDom<Domain,Dense>::add
  77. ( Rep& res, const Rep& u, const Rep& v ) const
  78. {
  79. GIVARO_ASSERT((dim(res) == dim(u)) && (dim(u) == dim(v)), "Bad size");
  80. AddOp<Domain> opcode (_domain);
  81. map( res, opcode, u, v);
  82. }
  83. template<class Domain>
  84. inline void VectorDom<Domain,Dense>::add
  85. ( Rep& res, const Rep& u, const Type_t& val ) const
  86. {
  87. GIVARO_ASSERT((dim(res) == dim(u)), "Bad size");
  88. Curried2<AddOp<Domain> > opcode(_domain, (Type_t&)val);
  89. map( res, opcode, u);
  90. }
  91. template<class Domain>
  92. inline void VectorDom<Domain,Dense>::add
  93. ( Rep& res, const Type_t& val, const Rep& v ) const
  94. {
  95. GIVARO_ASSERT((dim(res) == dim(v)), "Bad size");
  96. Curried1<AddOp<Domain> > opcode(_domain, val);
  97. map( res, opcode, v);
  98. }
  99. template<class Domain>
  100. inline void VectorDom<Domain,Dense>::add
  101. ( Rep& res, const VectorDom<Domain,Sparse>::Rep& u, const Rep& v ) const
  102. {
  103. GIVARO_ASSERT((dim(res) == dim(v)) && (dim(v) == u.dim()), "Bad size");
  104. // -- here assume : subdomains of res, u and v are equal
  105. assign(res, v);
  106. addin(res, u);
  107. }
  108. template<class Domain>
  109. inline void VectorDom<Domain,Dense>::add
  110. ( Rep& res, const Rep& v, const VectorDom<Domain,Sparse>::Rep& u ) const
  111. {
  112. GIVARO_ASSERT((dim(res) == dim(v)) && (dim(v) == u.dim()), "Bad size");
  113. // -- here assume : subdomains of res, u and v are equal
  114. assign(res, v);
  115. addin(res,u);
  116. }
  117. template<class Domain>
  118. inline void VectorDom<Domain,Dense>::addin
  119. ( Rep& res, const VectorDom<Domain,Sparse>::Rep& u ) const
  120. {
  121. GIVARO_ASSERT((dim(res) == u.dim()), "Bad size");
  122. // -- here assume : subdomains of res, u and v are equal
  123. typedef VectorDom<Domain,Sparse> VSparseDom;
  124. typedef VSparseDom::constIterator_t cIterator_t;
  125. typedef VSparseDom::IndiceIterator_t cIndiceIterator_t;
  126. cIterator_t u_curr = u.begin_data();
  127. cIterator_t u_end = u.end_data();
  128. cIndiceIterator_t u_indice = u.begin_indice();
  129. for ( ; u_curr != u_end; ++u_curr, ++u_indice )
  130. _domain.addin(res[*u_indice], *u_curr);
  131. }
  132. // --- sub
  133. template<class Domain>
  134. inline void VectorDom<Domain,Dense>::subin( Rep& res, const Rep& u ) const
  135. {
  136. GIVARO_ASSERT((dim(res) == dim(u)), "Bad size");
  137. SubOp<Domain> opcode ( _domain );
  138. map( res, opcode, u );
  139. }
  140. template<class Domain>
  141. inline void VectorDom<Domain,Dense>::sub( Rep& res, const Rep& u, const Rep& v ) const
  142. {
  143. GIVARO_ASSERT((dim(res) == dim(v)) && (dim(v) == dim(u)), "Bad size");
  144. SubOp<Domain> opcode ( _domain );
  145. map( res, opcode, u, v);
  146. }
  147. template<class Domain>
  148. inline void VectorDom<Domain,Dense>::sub( Rep& res, const Rep& u, const Type_t& val ) const
  149. {
  150. GIVARO_ASSERT((dim(res) == dim(u)), "Bad size");
  151. Curried2<SubOp<Domain> > opcode( _domain, (Type_t&)val);
  152. map( res, opcode, u);
  153. }
  154. template<class Domain>
  155. inline void VectorDom<Domain,Dense>::sub( Rep& res, const Type_t& val, const Rep& v ) const
  156. {
  157. GIVARO_ASSERT((dim(res) == dim(v)), "Bad size");
  158. Curried1<SubOp<Domain> > opcode( _domain, val);
  159. map( res, opcode, v);
  160. }
  161. template<class Domain>
  162. inline void VectorDom<Domain,Dense>::subin
  163. ( Rep& res, const VectorDom<Domain,Sparse>::Rep& u ) const
  164. {
  165. GIVARO_ASSERT((dim(res) == u.dim()), "Bad size");
  166. // -- here assume : subdomains of res, u and v are equal
  167. typedef VectorDom<Domain,Sparse> VSparseDom;
  168. typedef VSparseDom::constIterator_t cIterator_t;
  169. typedef VSparseDom::IndiceIterator_t cIndiceIterator_t;
  170. cIterator_t u_curr = u.begin_data();
  171. cIterator_t u_end = u.end_data();
  172. cIndiceIterator_t u_indice = u.begin_indice();
  173. for ( ; u_curr != u_end; ++u_curr, ++u_indice )
  174. _domain.subin(res[*u_indice], *u_curr);
  175. }
  176. template<class Domain>
  177. inline void VectorDom<Domain,Dense>::sub
  178. ( Rep& res, const VectorDom<Domain,Sparse>::Rep& u, const Rep& v ) const
  179. {
  180. GIVARO_ASSERT((dim(res) == dim(v)) && (dim(v) == u.dim()), "Bad size");
  181. // -- here assume : subdomains of res, u and v are equal
  182. neg(res, v);
  183. addin(res, u);
  184. }
  185. template<class Domain>
  186. inline void VectorDom<Domain,Dense>::sub
  187. ( Rep& res, const Rep& v, const VectorDom<Domain,Sparse>::Rep& u ) const
  188. {
  189. GIVARO_ASSERT((dim(res) == dim(v)) && (dim(v) == u.dim()), "Bad size");
  190. // -- here assume : subdomains of res, u and v are equal
  191. assign(res,v);
  192. subin(res,u);
  193. }
  194. // --- sub
  195. // --- neg
  196. template<class Domain>
  197. inline void VectorDom<Domain,Dense>::negin( Rep& res ) const
  198. {
  199. NegOp<Domain> opcode ( _domain );
  200. map( res, opcode, res );
  201. }
  202. template<class Domain>
  203. inline void VectorDom<Domain,Dense>::neg( Rep& res, const Rep& u ) const
  204. {
  205. GIVARO_ASSERT((dim(res) == dim(u)), "Bad size");
  206. NegOp<Domain> opcode ( _domain );
  207. map( res, opcode, u );
  208. }
  209. template<class Domain>
  210. void VectorDom<Domain,Dense>::axpy
  211. ( Rep& res, const Type_t& a, const Rep& x, const Rep& y ) const
  212. {
  213. GIVARO_ASSERT((dim(res) == dim(x)) && (dim(x) == dim(y)), "Bad size");
  214. size_t sz = dim(res);
  215. for ( size_t i=0; i<sz; ++i )
  216. _domain.axpy(res[i], a, x[i], y[i]);
  217. }
  218. template<class Domain>
  219. void VectorDom<Domain,Dense>::axpyin
  220. ( Rep& res, const Type_t& a, const Rep& x ) const
  221. {
  222. GIVARO_ASSERT((dim(res) == dim(x)), "Bad size");
  223. size_t sz = dim(res);
  224. for ( size_t i=0; i<sz; ++i )
  225. _domain.axpyin(res[i], a, x[i]);
  226. }
  227. template<class Domain>
  228. void VectorDom<Domain,Dense>::axmy
  229. ( Rep& res, const Type_t& a, const Rep& x, const Rep& y ) const
  230. {
  231. GIVARO_ASSERT((dim(res) == dim(x)) && (dim(x) == dim(y)), "Bad size");
  232. size_t sz = dim(res);
  233. for ( size_t i=0; i<sz; ++i )
  234. _domain.axmy(res[i], a, x[i], y[i]);
  235. }
  236. template<class Domain>
  237. void VectorDom<Domain,Dense>::axmyin
  238. ( Rep& res, const Type_t& a, const Rep& x ) const
  239. {
  240. GIVARO_ASSERT((dim(res) == dim(x)), "Bad size");
  241. size_t sz = dim(res);
  242. for ( size_t i=0; i<sz; ++i )
  243. _domain.axmyin(res[i], a, x[i]);
  244. }
  245. template<class Domain>
  246. void VectorDom<Domain,Dense>::axpy
  247. ( Rep& res, const Rep& a, const Rep& b, const Rep& y ) const
  248. {
  249. GIVARO_ASSERT((dim(res) == dim(a)) && (dim(a) == dim(b)) && (dim(b) == dim(y)), "Bad size");
  250. size_t sz = dim(res);
  251. for ( size_t i=0; i<sz; ++i )
  252. _domain.axpy(res[i], a[i], b[i], y[i]);
  253. }
  254. template<class Domain>
  255. void VectorDom<Domain,Dense>::axpyin
  256. ( Rep& res, const Rep& a, const Rep& b ) const
  257. {
  258. GIVARO_ASSERT((dim(res) == dim(a)) && (dim(a) == dim(b)), "Bad size");
  259. size_t sz = dim(res);
  260. for ( size_t i=0; i<sz; ++i )
  261. _domain.axpyin(res[i], a[i], b[i]);
  262. }
  263. template<class Domain>
  264. void VectorDom<Domain,Dense>::axmy
  265. ( Rep& res, const Rep& a, const Rep& b, const Rep& y ) const
  266. {
  267. GIVARO_ASSERT((dim(res) == dim(a)) && (dim(a) == dim(b)) && (dim(b) == dim(y)), "Bad size");
  268. size_t sz = dim(res);
  269. for ( size_t i=0; i<sz; ++i )
  270. _domain.axmy(res[i], a[i], b[i], y[i]);
  271. }
  272. template<class Domain>
  273. void VectorDom<Domain,Dense>::axmyin
  274. ( Rep& res, const Rep& a, const Rep& b ) const
  275. {
  276. GIVARO_ASSERT((dim(res) == dim(a)) && (dim(a) == dim(b)), "Bad size");
  277. size_t sz = dim(res);
  278. for ( size_t i=0; i<sz; ++i )
  279. _domain.axmyin(res[i], a[i], b[i]);
  280. }
  281. template<class Domain>
  282. inline void VectorDom<Domain,Dense>::dot
  283. ( Type_t& res, const Rep& op1, const Rep& op2) const
  284. {
  285. GIVARO_ASSERT((dim(op1) == dim(op2)), "Bad size");
  286. size_t sz = dim(op1);
  287. _domain.assign(res, _domain.zero);
  288. for ( size_t i=0; i<sz; ++i ) _domain.axpyin(res, op1[i], op2[i] );
  289. }
  290. // ==========================================================================
  291. //
  292. // -- Write the domain
  293. //
  294. template<class Domain>
  295. ostream& VectorDom<Domain, Dense>::write( ostream& o ) const
  296. {
  297. return _domain.write(o << "(") << ",Dense)";
  298. }
  299. // -- read the domain
  300. template<class Domain>
  301. istream& VectorDom<Domain, Dense>::read( istream& sin )
  302. {
  303. char ch;
  304. sin >> std::ws >> ch;
  305. if (ch != '(')
  306. GivError::throw_error(
  307. GivBadFormat("VectorDom<Domain,Dense>::read: syntax error no '('"));
  308. // -- read subdomain:
  309. _domain.read(sin);
  310. // -- read ,
  311. sin >> std::ws >> ch;
  312. if (ch != ',')
  313. GivError::throw_error(
  314. GivBadFormat("VectorDom<Domain,Dense>::read: syntax error no ','"));
  315. // -- read dense:
  316. char name[6];
  317. sin >> std::ws; sin.read(name, 5); name[5] = (char)0;
  318. if (strcmp(name, "Dense") !=0)
  319. GivError::throw_error(
  320. GivBadFormat("VectorDom<Domain,Dense>::read: syntax error no 'Dense'"));
  321. // -- read )
  322. sin >> std::ws >> ch;
  323. if (ch != ')')
  324. GivError::throw_error(
  325. GivBadFormat("VectorDom<Domain,Dense>::read: syntax error no ')'"));
  326. return sin;
  327. }
  328. // ==========================================================================
  329. //
  330. // -- Don't write the domain !, only write a rep / domain
  331. template<class Domain>
  332. ostream& VectorDom<Domain, Dense>::write(ostream& o, const Rep& V) const
  333. {
  334. if (dim(V) ==0) return o << "[ ]";
  335. o << '[';
  336. size_t i;
  337. for (i=0; i<dim(V)-1; i++) o << V[i] << ',';
  338. return o << V[i] << ']';
  339. }
  340. //
  341. // Read a vector given by a Maple-like format:
  342. // the grammar is :
  343. // s ---> '[' list_of_elt ']'
  344. // list_of_elt ---> elt
  345. // | list_of_elt ',' elt
  346. // The contraints are :
  347. // All lines are the same number of Elements.
  348. // The separators a those of the C lexical-convention, i.e.
  349. // ' ', '\n', '\t', '\f' .
  350. template<class Domain>
  351. istream& VectorDom<Domain,Dense>::read (istream& fin, Rep& A) const
  352. {
  353. char ch;
  354. // -- Skip the first "white":
  355. fin >> std::ws; fin.get(ch);
  356. if (ch != '[')
  357. GivError::throw_error(
  358. GivBadFormat("VectorDom<Domain,Dense>::read: syntax error no '['"));
  359. // -- Read the line and count the nb of elts
  360. int i = 0;
  361. Rep rep;
  362. rep.allocate(1);
  363. Type_t Tmp;
  364. fin >> Tmp;
  365. fin >> std::ws; fin.get(ch);
  366. rep[0] = Tmp;
  367. while (ch != ']')
  368. {
  369. if (ch != ',')
  370. GivError::throw_error(
  371. GivBadFormat("VectorDom<T,Dense>::read: syntax error no ','"));
  372. i++;
  373. fin >> std::ws >> Tmp;
  374. fin >> std::ws; fin.get(ch);
  375. // resize the vector :
  376. rep.reallocate(i+1);
  377. rep[i] = Tmp;
  378. }
  379. A.logcopy( rep );
  380. return fin;
  381. }
  382. } // Givaro
  383. #include "givaro/givvectdensespe.inl"