AbstractLinAlgPack: C++ Interfaces For Vectors, Matrices And Related Linear Algebra Objects  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
AbstractLinAlgPack_LinAlgOpPack.hpp
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Moocho: Multi-functional Object-Oriented arCHitecture for Optimization
5 // Copyright (2003) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Roscoe A. Bartlett (rabartl@sandia.gov)
38 //
39 // ***********************************************************************
40 // @HEADER
41 
42 #ifndef ABSTRACT_LIN_ALG_OP_PACK_H
43 #define ABSTRACT_LIN_ALG_OP_PACK_H
44 
45 #include "AbstractLinAlgPack_Types.hpp"
46 #include "AbstractLinAlgPack_MatrixOpNonsing.hpp"
47 #include "AbstractLinAlgPack_VectorMutable.hpp"
48 #include "AbstractLinAlgPack_VectorStdOps.hpp"
49 #include "AbstractLinAlgPack_AssertOp.hpp"
50 
51 namespace LinAlgOpPack {
52 
54 typedef AbstractLinAlgPack::size_type size_type;
56 typedef AbstractLinAlgPack::value_type value_type;
57 
70 
71 // Inject names of base linear algebra functions for DenseLinAlgPack.
72 // Note that this is neccesary in MS VC++ 5.0 because
73 // it does not perform name lookups properly but it
74 // is not adverse to the standard so it is a portable
75 // fix.
100 
117 
129 
134 template <class V>
135 void Vp_V(VectorMutable* v_lhs, const V& V_rhs);
136 
141 template <class V>
142 void assign(VectorMutable* v_lhs, const V& V_rhs);
143 
148 template <class V>
149 void V_StV(VectorMutable* v_lhs, value_type alpha, const V& V_rhs);
150 
155 template <class V>
156 void V_mV(VectorMutable* v_lhs, const V& V_rhs);
157 
159 
163 template <class V1, class V2>
164 void V_VpV(VectorMutable* v_lhs, const V1& V1_rhs1, const V2& V2_rhs2);
165 
170 template <class V1, class V2>
171 void V_VmV(VectorMutable* v_lhs, const V1& V1_rhs1, const V2& V2_rhs2);
172 
177 template <class V>
178 void V_StVpV(VectorMutable* v_lhs, value_type alpha, const V& V_rhs1
179  , const Vector& vs_rhs2);
180 
181 // end Level 1 BLAS for Vectors
183 
184 // //////////////////////////////////////////////////////////////////////////////
185 // ///////////////////////////////////////////////////////////////////////////////
197 
202 void Mp_M(MatrixOp* M_lhs, const MatrixOp& M_rhs, BLAS_Cpp::Transp trans_rhs);
203 
204 // end += operations
206 
211 void assign(MatrixOp* M_lhs, const MatrixOp& M_rhs, BLAS_Cpp::Transp trans_rhs);
212 
217 void M_StM(MatrixOp* M_lhs, value_type alpha, const MatrixOp& M_rhs
218  , BLAS_Cpp::Transp trans_rhs);
219 
224 void M_mM(MatrixOp* M_lhs, const MatrixOp& M_rhs, BLAS_Cpp::Transp trans_rhs) ;
225 
227 
231 void M_MpM(MatrixOp* M_lhs, const MatrixOp& M_rhs1, BLAS_Cpp::Transp trans_rhs1
232  , const MatrixOp& M_rhs2, BLAS_Cpp::Transp trans_rhs2);
233 
238 void M_MmM(MatrixOp* M_lhs, const MatrixOp& M_rhs1, BLAS_Cpp::Transp trans_rhs1
239  , const MatrixOp& M_rhs2, BLAS_Cpp::Transp trans_rhs2);
240 
245 void M_StMpM(MatrixOp* M_lhs, value_type alpha, const MatrixOp& M_rhs1, BLAS_Cpp::Transp trans_rhs1
246  , const MatrixOp& gms_rhs2, BLAS_Cpp::Transp trans_rhs2);
247 
248 // end Level 1 BLAS for Matrices
250 
251 // //////////////////////////////////////////////////////////////////////////////
252 // ///////////////////////////////////////////////////////////////////////
260 
265 template <class V>
266 void Vp_MtV(VectorMutable* v_lhs, const MatrixOp& M_rhs1, BLAS_Cpp::Transp trans_rhs1
267  , const V& V_rhs2);
268 
273 template <class V>
274 void V_StMtV(VectorMutable* v_lhs, value_type alpha, const MatrixOp& M_rhs1
275  , BLAS_Cpp::Transp trans_rhs1, const V& V_rhs2);
276 
281 template <class V>
282 void V_MtV(VectorMutable* v_lhs, const MatrixOp& M_rhs1, BLAS_Cpp::Transp trans_rhs1
283  , const V& V_rhs2);
284 
285 // end Level 2 BLAS
287 
288 // //////////////////////////////////////////////////////////////////////////////
289 // //////////////////////////////////////////////////////////////////////////////
297 
302 void Mp_MtM(MatrixOp* M_lhs, const MatrixOp& M_rhs1, BLAS_Cpp::Transp trans_rhs1
303  , const MatrixOp& M_rhs2, BLAS_Cpp::Transp trans_rhs2);
304 
309 void Mp_MtM(MatrixOp* M_lhs, const MatrixOp& M_rhs1, BLAS_Cpp::Transp trans_rhs1
310  , const MatrixOp& M_rhs2, BLAS_Cpp::Transp trans_rhs2, value_type beta);
311 
316 void M_StMtM(MatrixOp* M_lhs, value_type alpha, const MatrixOp& M_rhs1
317  , BLAS_Cpp::Transp trans_rhs1, const MatrixOp& M_rhs2, BLAS_Cpp::Transp trans_rhs2);
318 
323 void M_MtM(MatrixOp* M_lhs, const MatrixOp& M_rhs1
324  , BLAS_Cpp::Transp trans_rhs1, const MatrixOp& M_rhs2, BLAS_Cpp::Transp trans_rhs2);
325 
326 // end Level 3 BLAS
328 
329 // end Default Linear Algebra Implementations
331 
332 // ///////////////////////////////////////////////////////////////////////////
333 // ///////////////////////////////////////////////////////////////////////////
334 // Inline definitions
335 
336 // ///////////////////////////////////////////////////////////////////////////////
337 // Level 1 BLAS for Vectors
338 
339 // v_lhs += V_rhs.
340 template <class V>
341 inline
342 void Vp_V(VectorMutable* v_lhs, const V& V_rhs) {
343  Vp_StV(v_lhs,1.0,V_rhs);
344 }
345 
346 // v_lhs = - V_rhs.
347 template <class V>
348 inline
349 void V_mV(VectorMutable* v_lhs, const V& V_rhs) {
350  V_StV(v_lhs,-1.0,V_rhs);
351 }
352 
353 // ///////////////////////////////////////////////////////////////////////////////
354 // Level 1 BLAS for Matrices
355 
356 // M_lhs += op(M_rhs).
357 inline
358 void Mp_M(MatrixOp* M_lhs, const MatrixOp& M_rhs, BLAS_Cpp::Transp trans_rhs) {
359  Mp_StM(M_lhs,1.0,M_rhs,trans_rhs);
360 }
361 
362 // M_lhs = - op(M_rhs).
363 inline
364 void M_mM(MatrixOp* M_lhs, const MatrixOp& M_rhs, BLAS_Cpp::Transp trans_rhs) {
365  M_StM(M_lhs,-1.0,M_rhs,trans_rhs);
366 }
367 
368 // ///////////////////////////////////////////////////////////////////////
369 // Level 2 BLAS
370 
371 // v_lhs += op(M_rhs1) * V_rhs2.
372 template <class V>
373 inline
374 void Vp_MtV(VectorMutable* v_lhs, const MatrixOp& M_rhs1, BLAS_Cpp::Transp trans_rhs1
375  , const V& V_rhs2)
376 {
377  Vp_StMtV(v_lhs,1.0,M_rhs1,trans_rhs1,V_rhs2);
378 }
379 
380 // v_lhs = op(M_rhs1) * V_rhs2 + beta * v_lhs.
381 template <class V>
382 inline
383 void Vp_MtV(VectorMutable* v_lhs, const MatrixOp& M_rhs1, BLAS_Cpp::Transp trans_rhs1
384  , const V& V_rhs2, value_type beta)
385 {
386  Vp_StMtV(v_lhs,1.0,M_rhs1,trans_rhs1,V_rhs2,beta);
387 }
388 
389 // //////////////////////////////////////////////////////////////////////////////
390 // Level 3 BLAS
391 
392 // M_lhs += op(M_rhs1) * op(M_rhs2).
393 inline
394 void Mp_MtM(MatrixOp* M_lhs, const MatrixOp& M_rhs1, BLAS_Cpp::Transp trans_rhs1
395  , const MatrixOp& M_rhs2, BLAS_Cpp::Transp trans_rhs2)
396 {
397  Mp_StMtM(M_lhs,1.0,M_rhs1,trans_rhs1,M_rhs2,trans_rhs2);
398 }
399 
400 // M_lhs = op(M_rhs1) * op(M_rhs2) + beta * gms_rhs.
401 inline
402 void Mp_MtM(MatrixOp* M_lhs, const MatrixOp& M_rhs1, BLAS_Cpp::Transp trans_rhs1
403  , const MatrixOp& M_rhs2, BLAS_Cpp::Transp trans_rhs2, value_type beta)
404 {
405  Mp_StMtM(M_lhs,1.0,M_rhs1,trans_rhs1,M_rhs2,trans_rhs2,beta);
406 }
407 
408 // M_lhs = inv(op(M_rhs1)) * op(M_rhs2)
409 inline
410 void M_InvMtM(MatrixOp* M_lhs, const MatrixNonsing& M_rhs1, BLAS_Cpp::Transp trans_rhs1
411  , const MatrixOp& M_rhs2, BLAS_Cpp::Transp trans_rhs2 )
412 {
413  M_StInvMtM(M_lhs,1.0,M_rhs1,trans_rhs1,M_rhs2,trans_rhs2);
414 }
415 
416 } // end namespace LinAlgOpPack
417 
418 //
419 // Definitions
420 //
421 
422 namespace LinAlgOpPack {
423 
424 using BLAS_Cpp::rows;
425 using BLAS_Cpp::cols;
426 
427 // Inject assert functions
436 
437 // ///////////////////////////////////////////////////////////////////////////////
438 // Level 1 BLAS for Vectors
439 
440 // v_lhs = V_rhs.
441 template <class V>
442 void assign(VectorMutable* v_lhs, const V& V_rhs) {
443  Vp_V_assert_compatibility(v_lhs,V_rhs);
444  (*v_lhs) = 0.0;
445  Vp_V(v_lhs,V_rhs);
446 }
447 
448 // v_lhs = alpha * V_rhs.
449 template <class V>
450 void V_StV(VectorMutable* v_lhs, value_type alpha, const V& V_rhs) {
451  Vp_V_assert_compatibility(v_lhs,V_rhs);
452  (*v_lhs) = 0.0;
453  Vp_StV(v_lhs,alpha,V_rhs);
454 }
455 
456 // v_lhs = V1_rhs1 + V2_rhs2.
457 template <class V1, class V2>
458 void V_VpV(VectorMutable* v_lhs, const V1& V1_rhs1, const V2& V2_rhs2) {
459  VopV_assert_compatibility(V1_rhs1,V2_rhs2);
460  Vp_V_assert_compatibility(v_lhs,V1_rhs1);
461  (*v_lhs) = 0.0;
462  Vp_V(v_lhs,V1_rhs1);
463  Vp_V(v_lhs,V2_rhs2);
464 }
465 
466 // v_lhs = V_rhs1 - V_rhs2.
467 template <class V1, class V2>
468 void V_VmV(VectorMutable* v_lhs, const V1& V1_rhs1, const V2& V2_rhs2) {
469  VopV_assert_compatibility(V1_rhs1,V2_rhs2);
470  Vp_V_assert_compatibility(v_lhs,V1_rhs1);
471  (*v_lhs) = 0.0;
472  Vp_V(v_lhs,V1_rhs1);
473  Vp_StV(v_lhs,-1.0,V2_rhs2);
474 }
475 
476 // v_lhs = alpha * V_rhs1 + v_rhs2.
477 template <class V>
478 void V_StVpV(VectorMutable* v_lhs, value_type alpha, const V& V_rhs1
479  , const Vector& v_rhs2)
480 {
481  VopV_assert_compatibility(V_rhs1,v_rhs2);
482  (*v_lhs) = v_rhs2;
483  Vp_StV(v_lhs,alpha,V_rhs1);
484 }
485 
486 // ////////////////////////////////////////////////////////////////
487 // Level 2 BLAS
488 
489 // v_lhs = alpha * op(M_rhs1) * V_rhs2.
490 template <class V>
491 void V_StMtV(VectorMutable* v_lhs, value_type alpha, const MatrixOp& M_rhs1
492  , BLAS_Cpp::Transp trans_rhs1, const V& V_rhs2)
493 {
494  Vp_MtV_assert_compatibility(v_lhs,M_rhs1,trans_rhs1,V_rhs2);
495  Vp_StMtV(v_lhs,alpha,M_rhs1,trans_rhs1,V_rhs2,0.0);
496 }
497 
498 // v_lhs = op(M_rhs1) * V_rhs2.
499 template <class V>
500 void V_MtV(VectorMutable* v_lhs, const MatrixOp& M_rhs1, BLAS_Cpp::Transp trans_rhs1
501  , const V& V_rhs2)
502 {
503  Vp_MtV_assert_compatibility(v_lhs,M_rhs1,trans_rhs1,V_rhs2);
504  Vp_StMtV(v_lhs,1.0,M_rhs1,trans_rhs1,V_rhs2,0.0);
505 }
506 
507 } // end namespace LinAlgOpPack
508 
509 #endif // ABSTRACT_LIN_ALG_OP_PACK_H
void V_StVpV(VectorMutable *v_lhs, value_type alpha, const V &V_rhs1, const Vector &vs_rhs2)
v_lhs = alpha * V_rhs1 + vs_rhs2.
Abstract interface for immutable, finite dimensional, coordinate vectors {abstract}.
void Mp_StMtM(MatrixOp *mwo_lhs, value_type alpha, const MatrixOp &mwo_rhs1, BLAS_Cpp::Transp trans_rhs1, const MatrixOp &mwo_rhs2, BLAS_Cpp::Transp trans_rhs2, value_type beta=1.0)
mwo_lhs = alpha * op(mwo_rhs1) * op(mwo_rhs2) + beta * mwo_lhs (right) (xGEMM).
void Vt_S(VectorMutable *v_lhs, const value_type &alpha)
v_lhs *= alpha
void Vp_StV(VectorMutable *v_lhs, const value_type &alpha, const Vector &v_rhs)
v_lhs = alpha * v_rhs + v_lhs
size_type rows(size_type rows, size_type cols, BLAS_Cpp::Transp _trans)
void M_mM(MatrixOp *M_lhs, const MatrixOp &M_rhs, BLAS_Cpp::Transp trans_rhs)
gm_lhs = - op(M_rhs).
void V_StV(VectorMutable *v_lhs, value_type alpha, const V &V_rhs)
v_lhs = alpha * V_rhs.
void MopM_assert_compatibility(const MatrixOp &m_rhs1, BLAS_Cpp::Transp trans_rhs1, const MatrixOp &m_rhs2, BLAS_Cpp::Transp trans_rhs2)
op(m_rhs1) op op(m_rhs2)
void MtV_assert_compatibility(const MatrixOp &m_rhs1, BLAS_Cpp::Transp trans_rhs1, const Vector &v_rhs2)
op(m_rhs1) * v_rhs2
void Vp_V_assert_compatibility(VectorMutable *v_lhs, const Vector &v_rhs)
v_lhs += op v_rhs
void M_MmM(MatrixOp *M_lhs, const MatrixOp &M_rhs1, BLAS_Cpp::Transp trans_rhs1, const MatrixOp &M_rhs2, BLAS_Cpp::Transp trans_rhs2)
M_lhs = op(M_rhs1) - op(M_rhs2).
void assign(MatrixOp *M_lhs, const MatrixOp &M_rhs, BLAS_Cpp::Transp trans_rhs)
M_lhs = op(M_rhs).
void V_mV(VectorMutable *v_lhs, const V &V_rhs)
v_lhs = - V_rhs.
void MtM_assert_compatibility(const MatrixOp &m_rhs1, BLAS_Cpp::Transp trans_rhs1, const MatrixOp &m_rhs2, BLAS_Cpp::Transp trans_rhs2)
op(m_lhs) += op(m_rhs1)
Abstract interface for objects that represent a space for mutable coordinate vectors.
void M_MpM(MatrixOp *M_lhs, const MatrixOp &M_rhs1, BLAS_Cpp::Transp trans_rhs1, const MatrixOp &M_rhs2, BLAS_Cpp::Transp trans_rhs2)
void M_StM(MatrixOp *M_lhs, value_type alpha, const MatrixOp &M_rhs, BLAS_Cpp::Transp trans_rhs)
gm_lhs = alpha * M_rhs.
void V_VpV(VectorMutable *v_lhs, const V1 &V1_rhs1, const V2 &V2_rhs2)
size_t size_type
void syrk(const MatrixOp &mwo_rhs, BLAS_Cpp::Transp M_trans, value_type alpha, value_type beta, MatrixSymOp *sym_lhs)
Perform a rank-k update of a symmetric matrix of the form:
void Mp_MtM_assert_compatibility(MatrixOp *m_lhs, BLAS_Cpp::Transp trans_lhs, const MatrixOp &m_rhs1, BLAS_Cpp::Transp trans_rhs1, const MatrixOp &m_rhs2, BLAS_Cpp::Transp trans_rhs2)
op(m_lhs) += op(m_rhs1) * op(m_rhs2)
void V_InvMtV(VectorMutable *v_lhs, const MatrixNonsing &M_rhs1, BLAS_Cpp::Transp trans_rhs1, const Vector &v_rhs2)
v_lhs = inv(op(M_rhs1)) * v_rhs2
void Mp_StM(MatrixOp *mwo_lhs, value_type alpha, const MatrixOp &M_rhs, BLAS_Cpp::Transp trans_rhs)
void Vp_StMtV(VectorMutable *v_lhs, value_type alpha, const MatrixOp &M_rhs1, BLAS_Cpp::Transp trans_rhs1, const Vector &v_rhs2, value_type beta=1.0)
v_lhs = alpha * op(M_rhs1) * v_rhs2 + beta * v_lhs (BLAS xGEMV)
void Mp_M(MatrixOp *M_lhs, const MatrixOp &M_rhs, BLAS_Cpp::Transp trans_rhs)
M_lhs += op(M_rhs).
Base class for all matrices that support basic matrix operations.
void Mt_S(MatrixOp *mwo_lhs, value_type alpha)
value_type dot(const Vector &v_rhs1, const Vector &v_rhs2)
result = v_rhs1' * v_rhs2
Abstract interface for mutable coordinate vectors {abstract}.
void Vp_S(VectorMutable *v_lhs, const value_type &alpha)
v_lhs += alpha
value_type sum(const Vector &v_rhs)
result = sum( v_rhs(i), i = 1,,,dim )
Abstract base class for all nonsingular polymorphic matrices that can be used to compute matrix-vecto...
void Vp_MtV_assert_compatibility(VectorMutable *v_lhs, const MatrixOp &m_rhs1, BLAS_Cpp::Transp trans_rhs1, const Vector &v_rhs2)
v_lhs += op(m_rhs1) * v_rhs2
void M_StInvMtM(MatrixOp *m_lhs, value_type alpha, const MatrixNonsing &M_rhs1, BLAS_Cpp::Transp trans_rhs1, const MatrixOp &mwo_rhs2, BLAS_Cpp::Transp trans_rhs2)
m_lhs = alpha * inv(op(mwo_rhs1)) * op(mwo_rhs2) (right)
void M_StMpM(MatrixOp *M_lhs, value_type alpha, const MatrixOp &M_rhs1, BLAS_Cpp::Transp trans_rhs1, const MatrixOp &gms_rhs2, BLAS_Cpp::Transp trans_rhs2)
M_lhs = alpha * op(M_rhs1) + op(gms_rhs2).
Abstract base class for all nonsingular polymorphic matrices that can solve for linear system with bu...
Transp
size_type cols(size_type rows, size_type cols, BLAS_Cpp::Transp _trans)
void V_VmV(VectorMutable *v_lhs, const V1 &V1_rhs1, const V2 &V2_rhs2)
v_lhs = V_rhs1 - V_rhs2.
void Mp_M_assert_compatibility(MatrixOp *m_lhs, BLAS_Cpp::Transp trans_lhs, const MatrixOp &m_rhs, BLAS_Cpp::Transp trans_rhs)
op(m_lhs) += op op(m_rhs)
void Vp_V(VectorMutable *v_lhs, const V &V_rhs)
v_lhs += V_rhs.
void VopV_assert_compatibility(const Vector &v_rhs1, const Vector &v_rhs2)
v_rhs1 op v_rhs2