Warning: this documentation for the development version is under construction.
00001 // Copyright (C) 2001-2009 Vivien Mallet 00002 // 00003 // This file is part of the linear-algebra library Seldon, 00004 // http://seldon.sourceforge.net/. 00005 // 00006 // Seldon is free software; you can redistribute it and/or modify it under the 00007 // terms of the GNU Lesser General Public License as published by the Free 00008 // Software Foundation; either version 2.1 of the License, or (at your option) 00009 // any later version. 00010 // 00011 // Seldon is distributed in the hope that it will be useful, but WITHOUT ANY 00012 // WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS 00013 // FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for 00014 // more details. 00015 // 00016 // You should have received a copy of the GNU Lesser General Public License 00017 // along with Seldon. If not, see http://www.gnu.org/licenses/. 00018 00019 00020 #ifndef SELDON_FILE_BLAS_1_CXX 00021 00022 00023 #include "Blas_1.hxx" 00024 00025 00026 namespace Seldon 00027 { 00028 00029 00031 // GenRot // 00032 00033 00034 void GenRot(float& a, float& b, float& c, float& d) 00035 { 00036 cblas_srotg(&a, &b, &c, &d); 00037 } 00038 00039 00040 void GenRot(double& a, double& b, double& c, double& d) 00041 { 00042 cblas_drotg(&a, &b, &c, &d); 00043 } 00044 00045 00046 // GenRot // 00048 00049 00050 00052 // GenModifRot // 00053 00054 00055 void GenModifRot(float& d1, float& d2, 00056 float& x1, const float& y1, 00057 float* param) 00058 { 00059 cblas_srotmg(&d1, &d2, &x1, y1, param); 00060 } 00061 00062 00063 void GenModifRot(double& d1, double& d2, 00064 double& x1, const double& y1, 00065 double* param) 00066 { 00067 cblas_drotmg(&d1, &d2, &x1, y1, param); 00068 } 00069 00070 00071 // GenModifRot // 00073 00074 00075 00077 // ApplyRot // 00078 00079 00080 template <class Allocator> 00081 void ApplyRot(Vector<float, VectFull, Allocator>& X, 00082 Vector<float, VectFull, Allocator>& Y, 00083 const float c, const float s) 00084 { 00085 cblas_srot(X.GetLength(), X.GetData(), 1, 00086 Y.GetData(), 1, c, s); 00087 } 00088 00089 00090 template <class Allocator> 00091 void ApplyRot(Vector<double, VectFull, Allocator>& X, 00092 Vector<double, VectFull, Allocator>& Y, 00093 const double c, const double s) 00094 { 00095 cblas_drot(X.GetLength(), X.GetData(), 1, 00096 Y.GetData(), 1, c, s); 00097 } 00098 00099 00100 // ApplyRot // 00102 00103 00104 00106 // ApplyModifRot // 00107 00108 00109 template <class Allocator> 00110 void ApplyModifRot(Vector<float, VectFull, Allocator>& X, 00111 Vector<float, VectFull, Allocator>& Y, 00112 const float* param) 00113 { 00114 cblas_srotm(X.GetLength(), X.GetData(), 1, 00115 Y.GetData(), 1, param); 00116 } 00117 00118 00119 template <class Allocator> 00120 void ApplyModifRot(Vector<double, VectFull, Allocator>& X, 00121 Vector<double, VectFull, Allocator>& Y, 00122 const double* param) 00123 { 00124 cblas_drotm(X.GetLength(), X.GetData(), 1, 00125 Y.GetData(), 1, param); 00126 } 00127 00128 00129 // ApplyModifRot // 00131 00132 00133 00135 // Swap // 00136 00137 00138 template <class Allocator> 00139 void Swap(Vector<float, VectFull, Allocator>& X, 00140 Vector<float, VectFull, Allocator>& Y) 00141 { 00142 00143 #ifdef SELDON_CHECK_DIMENSIONS 00144 CheckDim(X, Y, "Swap(X, Y)", "X <-> Y"); 00145 #endif 00146 00147 cblas_sswap(X.GetLength(), X.GetData(), 1, 00148 Y.GetData(), 1); 00149 } 00150 00151 00152 template <class Allocator> 00153 void Swap(Vector<double, VectFull, Allocator>& X, 00154 Vector<double, VectFull, Allocator>& Y) 00155 { 00156 00157 #ifdef SELDON_CHECK_DIMENSIONS 00158 CheckDim(X, Y, "Swap(X, Y)", "X <-> Y"); 00159 #endif 00160 00161 cblas_dswap(X.GetLength(), X.GetData(), 1, 00162 Y.GetData(), 1); 00163 } 00164 00165 00166 template <class Allocator> 00167 void Swap(Vector<complex<float>, VectFull, Allocator>& X, 00168 Vector<complex<float>, VectFull, Allocator>& Y) 00169 { 00170 00171 #ifdef SELDON_CHECK_DIMENSIONS 00172 CheckDim(X, Y, "Swap(X, Y)", "X <-> Y"); 00173 #endif 00174 00175 cblas_cswap(X.GetLength(), reinterpret_cast<void*>(X.GetData()), 1, 00176 reinterpret_cast<void*>(Y.GetData()), 1); 00177 } 00178 00179 00180 template <class Allocator> 00181 void Swap(Vector<complex<double>, VectFull, Allocator>& X, 00182 Vector<complex<double>, VectFull, Allocator>& Y) 00183 { 00184 00185 #ifdef SELDON_CHECK_DIMENSIONS 00186 CheckDim(X, Y, "Swap(X, Y)", "X <-> Y"); 00187 #endif 00188 00189 cblas_zswap(X.GetLength(), reinterpret_cast<void*>(X.GetData()), 1, 00190 reinterpret_cast<void*>(Y.GetData()), 1); 00191 } 00192 00193 00194 // Swap // 00196 00197 00198 00200 // Mlt // 00201 00202 00203 template <class Allocator> 00204 void Mlt(const float alpha, 00205 Vector<float, VectFull, Allocator>& X) 00206 { 00207 cblas_sscal(X.GetLength(), alpha, X.GetData(), 1); 00208 } 00209 00210 00211 template <class Allocator> 00212 void Mlt(const double alpha, 00213 Vector<double, VectFull, Allocator>& X) 00214 { 00215 cblas_dscal(X.GetLength(), alpha, X.GetData(), 1); 00216 } 00217 00218 00219 template <class Allocator> 00220 void Mlt(const float alpha, 00221 Vector<complex<float>, VectFull, Allocator>& X) 00222 { 00223 cblas_csscal(X.GetLength(), alpha, 00224 reinterpret_cast<void*>(X.GetData()), 1); 00225 } 00226 00227 00228 template <class Allocator> 00229 void Mlt(const double alpha, 00230 Vector<complex<double>, VectFull, Allocator>& X) 00231 { 00232 cblas_zdscal(X.GetLength(), alpha, 00233 reinterpret_cast<void*>(X.GetData()), 1); 00234 } 00235 00236 00237 template <class Allocator> 00238 void Mlt(const complex<float> alpha, 00239 Vector<complex<float>, VectFull, Allocator>& X) 00240 { 00241 cblas_cscal(X.GetLength(), 00242 reinterpret_cast<const void*>(&alpha), 00243 reinterpret_cast<void*>(X.GetData()), 1); 00244 } 00245 00246 00247 template <class Allocator> 00248 void Mlt(const complex<double> alpha, 00249 Vector<complex<double>, VectFull, Allocator>& X) 00250 { 00251 cblas_zscal(X.GetLength(), 00252 reinterpret_cast<const void*>(&alpha), 00253 reinterpret_cast<void*>(X.GetData()), 1); 00254 } 00255 00256 00257 // Mlt // 00259 00260 00261 00263 // Copy // 00264 00265 00266 template <class Allocator0, class Allocator1> 00267 void Copy(const Vector<float, VectFull, Allocator0>& X, 00268 Vector<float, VectFull, Allocator1>& Y) 00269 { 00270 00271 #ifdef SELDON_CHECK_DIMENSIONS 00272 CheckDim(X, Y, "Copy(X, Y)", "X -> Y"); 00273 #endif 00274 00275 cblas_scopy(Y.GetLength(), 00276 reinterpret_cast<const float*>(X.GetData()), 1, 00277 reinterpret_cast<float*>(Y.GetData()), 1); 00278 } 00279 00280 00281 template <class Allocator0, class Allocator1> 00282 void Copy(const Vector<double, VectFull, Allocator0>& X, 00283 Vector<double, VectFull, Allocator1>& Y) 00284 { 00285 00286 #ifdef SELDON_CHECK_DIMENSIONS 00287 CheckDim(X, Y, "Copy(X, Y)", "X -> Y"); 00288 #endif 00289 00290 cblas_dcopy(Y.GetLength(), 00291 reinterpret_cast<const double*>(X.GetData()), 1, 00292 reinterpret_cast<double*>(Y.GetData()), 1); 00293 } 00294 00295 00296 template <class Allocator0, class Allocator1> 00297 void Copy(const Vector<complex<float>, VectFull, Allocator0>& X, 00298 Vector<complex<float>, VectFull, Allocator1>& Y) 00299 { 00300 00301 #ifdef SELDON_CHECK_DIMENSIONS 00302 CheckDim(X, Y, "Copy(X, Y)", "X -> Y"); 00303 #endif 00304 00305 cblas_ccopy(Y.GetLength(), 00306 reinterpret_cast<const void*>(X.GetData()), 1, 00307 reinterpret_cast<void*>(Y.GetData()), 1); 00308 } 00309 00310 00311 template <class Allocator0, class Allocator1> 00312 void Copy(const Vector<complex<double>, VectFull, Allocator0>& X, 00313 Vector<complex<double>, VectFull, Allocator1>& Y) 00314 { 00315 00316 #ifdef SELDON_CHECK_DIMENSIONS 00317 CheckDim(X, Y, "Copy(X, Y)", "X -> Y"); 00318 #endif 00319 00320 cblas_zcopy(Y.GetLength(), 00321 reinterpret_cast<const void*>(X.GetData()), 1, 00322 reinterpret_cast<void*>(Y.GetData()), 1); 00323 } 00324 00325 00326 // Copy // 00328 00329 00330 00332 // Add // 00333 00334 00335 template <class Allocator0, class Allocator1> 00336 void Add(const float alpha, 00337 const Vector<float, VectFull, Allocator0>& X, 00338 Vector<float, VectFull, Allocator1>& Y) 00339 { 00340 00341 #ifdef SELDON_CHECK_DIMENSIONS 00342 CheckDim(X, Y, "Add(alpha, X, Y)"); 00343 #endif 00344 00345 cblas_saxpy(Y.GetLength(), 00346 alpha, 00347 reinterpret_cast<const float*>(X.GetData()), 1, 00348 reinterpret_cast<float*>(Y.GetData()), 1); 00349 } 00350 00351 00352 template <class Allocator0, class Allocator1> 00353 void Add(const double alpha, 00354 const Vector<double, VectFull, Allocator0>& X, 00355 Vector<double, VectFull, Allocator1>& Y) 00356 { 00357 00358 #ifdef SELDON_CHECK_DIMENSIONS 00359 CheckDim(X, Y, "Add(alpha, X, Y)"); 00360 #endif 00361 00362 cblas_daxpy(Y.GetLength(), 00363 alpha, 00364 reinterpret_cast<const double*>(X.GetData()), 1, 00365 reinterpret_cast<double*>(Y.GetData()), 1); 00366 } 00367 00368 00369 template <class Allocator0, class Allocator1> 00370 void Add(const complex<float> alpha, 00371 const Vector<complex<float>, VectFull, Allocator0>& X, 00372 Vector<complex<float>, VectFull, Allocator1>& Y) 00373 { 00374 00375 #ifdef SELDON_CHECK_DIMENSIONS 00376 CheckDim(X, Y, "Add(alpha, X, Y)"); 00377 #endif 00378 00379 cblas_caxpy(Y.GetLength(), 00380 reinterpret_cast<const void*>(&alpha), 00381 reinterpret_cast<const void*>(X.GetData()), 1, 00382 reinterpret_cast<float*>(Y.GetData()), 1); 00383 } 00384 00385 00386 template <class Allocator0, class Allocator1> 00387 void Add(const complex<double> alpha, 00388 const Vector<complex<double>, VectFull, Allocator0>& X, 00389 Vector<complex<double>, VectFull, Allocator1>& Y) 00390 { 00391 00392 #ifdef SELDON_CHECK_DIMENSIONS 00393 CheckDim(X, Y, "Add(alpha, X, Y)"); 00394 #endif 00395 00396 cblas_zaxpy(Y.GetLength(), 00397 reinterpret_cast<const void*>(&alpha), 00398 reinterpret_cast<const void*>(X.GetData()), 1, 00399 reinterpret_cast<double*>(Y.GetData()), 1); 00400 } 00401 00402 00403 // Add // 00405 00406 00407 00409 // DotProd // 00410 00411 00412 template <class Allocator0, class Allocator1> 00413 float DotProd(const Vector<float, VectFull, Allocator0>& X, 00414 const Vector<float, VectFull, Allocator1>& Y) 00415 { 00416 00417 #ifdef SELDON_CHECK_DIMENSIONS 00418 CheckDim(X, Y, "DotProd(X, Y)", "dot(X, Y)"); 00419 #endif 00420 00421 return cblas_sdot(Y.GetLength(), 00422 reinterpret_cast<const float*>(X.GetData()), 1, 00423 reinterpret_cast<const float*>(Y.GetData()), 1); 00424 } 00425 00426 00427 template <class Allocator0, class Allocator1> 00428 double DotProd(const Vector<double, VectFull, Allocator0>& X, 00429 const Vector<double, VectFull, Allocator1>& Y) 00430 { 00431 00432 #ifdef SELDON_CHECK_DIMENSIONS 00433 CheckDim(X, Y, "DotProd(X, Y)", "dot(X, Y)"); 00434 #endif 00435 00436 return cblas_ddot(Y.GetLength(), 00437 reinterpret_cast<const double*>(X.GetData()), 1, 00438 reinterpret_cast<const double*>(Y.GetData()), 1); 00439 } 00440 00441 00442 template <class Allocator0, class Allocator1> 00443 complex<float> 00444 DotProd(const Vector<complex<float>, VectFull, Allocator0>& X, 00445 const Vector<complex<float>, VectFull, Allocator1>& Y) 00446 { 00447 00448 #ifdef SELDON_CHECK_DIMENSIONS 00449 CheckDim(X, Y, "DotProd(X, Y)", "dot(X, Y)"); 00450 #endif 00451 00452 complex<float> dotu; 00453 cblas_cdotu_sub(Y.GetLength(), 00454 reinterpret_cast<const void*>(X.GetData()), 1, 00455 reinterpret_cast<const void*>(Y.GetData()), 1, 00456 reinterpret_cast<void*>(&dotu)); 00457 return dotu; 00458 } 00459 00460 00461 template <class Allocator0, class Allocator1> 00462 complex<double> 00463 DotProd(const Vector<complex<double>, VectFull, Allocator0>& X, 00464 const Vector<complex<double>, VectFull, Allocator1>& Y) 00465 { 00466 00467 #ifdef SELDON_CHECK_DIMENSIONS 00468 CheckDim(X, Y, "DotProd(X, Y)", "dot(X, Y)"); 00469 #endif 00470 00471 complex<double> dotu; 00472 cblas_zdotu_sub(Y.GetLength(), 00473 reinterpret_cast<const void*>(X.GetData()), 1, 00474 reinterpret_cast<const void*>(Y.GetData()), 1, 00475 reinterpret_cast<void*>(&dotu)); 00476 return dotu; 00477 } 00478 00479 00480 // DotProd // 00482 00483 00484 00486 // SCALEDDOTPROD // 00487 00488 00489 template <class Allocator0, class Allocator1> 00490 float ScaledDotProd(const float alpha, 00491 const Vector<float, VectFull, Allocator0>& X, 00492 const Vector<float, VectFull, Allocator1>& Y) 00493 { 00494 00495 #ifdef SELDON_CHECK_DIMENSIONS 00496 CheckDim(X, Y, "ScaledDotProd(X, Y)", "dot(X, Y)"); 00497 #endif 00498 00499 return cblas_sdsdot(Y.GetLength(), alpha, 00500 reinterpret_cast<const float*>(X.GetData()), 1, 00501 reinterpret_cast<const float*>(Y.GetData()), 1); 00502 } 00503 00504 00505 // SCALEDDOTPROD // 00507 00508 00509 00511 // DotProjConj // 00512 00513 00514 template <class Allocator0, class Allocator1> 00515 complex<float> 00516 DotProdConj(const Vector<complex<float>, VectFull, Allocator0>& X, 00517 const Vector<complex<float>, VectFull, Allocator1>& Y) 00518 { 00519 00520 #ifdef SELDON_CHECK_DIMENSIONS 00521 CheckDim(X, Y, "DotProdConj(X, Y)", "dot(X, Y)"); 00522 #endif 00523 00524 complex<float> dotc; 00525 cblas_cdotc_sub(Y.GetLength(), 00526 reinterpret_cast<const void*>(X.GetData()), 1, 00527 reinterpret_cast<const void*>(Y.GetData()), 1, 00528 reinterpret_cast<void*>(&dotc)); 00529 return dotc; 00530 } 00531 00532 00533 template <class Allocator0, class Allocator1> 00534 complex<double> 00535 DotProdConj(const Vector<complex<double>, VectFull, Allocator0>& X, 00536 const Vector<complex<double>, VectFull, Allocator1>& Y) 00537 { 00538 00539 #ifdef SELDON_CHECK_DIMENSIONS 00540 CheckDim(X, Y, "DotProdConj(X, Y)", "dot(X, Y)"); 00541 #endif 00542 00543 complex<double> dotc; 00544 cblas_zdotc_sub(Y.GetLength(), 00545 reinterpret_cast<const void*>(X.GetData()), 1, 00546 reinterpret_cast<const void*>(Y.GetData()), 1, 00547 reinterpret_cast<void*>(&dotc)); 00548 return dotc; 00549 } 00550 00551 00552 // DotProdConj // 00554 00555 00556 00558 // Norm1 // 00559 00560 00561 template <class Allocator> 00562 float Norm1(const Vector<float, VectFull, Allocator>& X) 00563 { 00564 return cblas_sasum(X.GetLength(), 00565 reinterpret_cast<const float*>(X.GetData()), 1); 00566 } 00567 00568 00569 template <class Allocator> 00570 double Norm1(const Vector<double, VectFull, Allocator>& X) 00571 { 00572 return cblas_dasum(X.GetLength(), 00573 reinterpret_cast<const double*>(X.GetData()), 1); 00574 } 00575 00576 00577 template <class Allocator> 00578 float Norm1(const Vector<complex<float>, VectFull, Allocator>& X) 00579 { 00580 return cblas_scasum(X.GetLength(), 00581 reinterpret_cast<const void*>(X.GetData()), 1); 00582 } 00583 00584 00585 template <class Allocator> 00586 double Norm1(const Vector<complex<double>, VectFull, Allocator>& X) 00587 { 00588 return cblas_dzasum(X.GetLength(), 00589 reinterpret_cast<const void*>(X.GetData()), 1); 00590 } 00591 00592 00593 // Norm1 // 00595 00596 00597 00599 // Norm2 // 00600 00601 00602 template <class Allocator> 00603 float Norm2(const Vector<float, VectFull, Allocator>& X) 00604 { 00605 return cblas_snrm2(X.GetLength(), 00606 reinterpret_cast<const float*>(X.GetData()), 1); 00607 } 00608 00609 00610 template <class Allocator> 00611 double Norm2(const Vector<double, VectFull, Allocator>& X) 00612 { 00613 return cblas_dnrm2(X.GetLength(), 00614 reinterpret_cast<const double*>(X.GetData()), 1); 00615 } 00616 00617 00618 template <class Allocator> 00619 float Norm2(const Vector<complex<float>, VectFull, Allocator>& X) 00620 { 00621 return cblas_scnrm2(X.GetLength(), 00622 reinterpret_cast<const void*>(X.GetData()), 1); 00623 } 00624 00625 00626 template <class Allocator> 00627 double Norm2(const Vector<complex<double>, VectFull, Allocator>& X) 00628 { 00629 return cblas_dznrm2(X.GetLength(), 00630 reinterpret_cast<const void*>(X.GetData()), 1); 00631 } 00632 00633 00634 // Norm2 // 00636 00637 00638 00640 // GetMaxAbsIndex // 00641 00642 00643 template <class Allocator> 00644 size_t GetMaxAbsIndex(const Vector<float, VectFull, Allocator>& X) 00645 { 00646 return cblas_isamax(X.GetLength(), 00647 reinterpret_cast<const float*>(X.GetData()), 1); 00648 } 00649 00650 00651 template <class Allocator> 00652 size_t GetMaxAbsIndex(const Vector<double, VectFull, Allocator>& X) 00653 { 00654 return cblas_idamax(X.GetLength(), 00655 reinterpret_cast<const double*>(X.GetData()), 1); 00656 } 00657 00658 00659 template <class Allocator> 00660 size_t GetMaxAbsIndex(const Vector<complex<float>, VectFull, Allocator>& X) 00661 { 00662 return cblas_icamax(X.GetLength(), 00663 reinterpret_cast<const void*>(X.GetData()), 1); 00664 } 00665 00666 00667 template <class Allocator> 00668 size_t 00669 GetMaxAbsIndex(const Vector<complex<double>, VectFull, Allocator>& X) 00670 { 00671 return cblas_izamax(X.GetLength(), 00672 reinterpret_cast<const void*>(X.GetData()), 1); 00673 } 00674 00675 00676 // GetMaxAbsIndex // 00678 00679 00680 } // namespace Seldon. 00681 00682 #define SELDON_FILE_BLAS_1_CXX 00683 #endif