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_MatrixNonsingSerial.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 SLAP_MATRIX_NONSINGULAR_SERIAL_H
43 #define SLAP_MATRIX_NONSINGULAR_SERIAL_H
44 
45 #include "AbstractLinAlgPack_Types.hpp"
46 #include "AbstractLinAlgPack_MatrixNonsing.hpp"
47 
48 namespace AbstractLinAlgPack {
49 
60  : virtual public AbstractLinAlgPack::MatrixNonsing // doxygen needs full name
61 {
62 public:
63 
66 
68  virtual void V_InvMtV(
69  DVector* v_lhs, BLAS_Cpp::Transp trans_rhs1
70  ,const DVectorSlice& vs_rhs2) const;
72  virtual void V_InvMtV(
73  DVectorSlice* vs_lhs, BLAS_Cpp::Transp trans_rhs1
74  ,const DVectorSlice& vs_rhs2) const = 0;
76  virtual void V_InvMtV(
77  DVector* v_lhs, BLAS_Cpp::Transp trans_rhs1
78  ,const SpVectorSlice& sv_rhs2) const;
80  virtual void V_InvMtV(
81  DVectorSlice* vs_lhs, BLAS_Cpp::Transp trans_rhs1
82  ,const SpVectorSlice& sv_rhs2) const;
84  virtual value_type transVtInvMtV(
85  const DVectorSlice& vs_rhs1, BLAS_Cpp::Transp trans_rhs2, const DVectorSlice& vs_rhs3) const;
87  virtual value_type transVtInvMtV(
88  const SpVectorSlice& sv_rhs1, BLAS_Cpp::Transp trans_rhs2, const SpVectorSlice& sv_rhs3) const;
89 
90  // end Level-2 BLAS
92 
95 
97  virtual void M_StInvMtM(
98  DMatrix* gm_lhs, value_type alpha
99  ,BLAS_Cpp::Transp trans_rhs1
100  ,const DMatrixSlice& gms_rhs2, BLAS_Cpp::Transp trans_rhs2 ) const;
102  virtual void M_StInvMtM(
103  DMatrixSlice* gms_lhs, value_type alpha
104  ,BLAS_Cpp::Transp trans_rhs1
105  ,const DMatrixSlice& gms_rhs2, BLAS_Cpp::Transp trans_rhs2 ) const;
107  virtual void M_StMtInvM(
108  DMatrix* gm_lhs, value_type alpha
109  ,const DMatrixSlice& gms_rhs1, BLAS_Cpp::Transp trans_rhs1
110  ,BLAS_Cpp::Transp trans_rhs2 ) const;
112  virtual void M_StMtInvM(
113  DMatrixSlice* gms_lhs, value_type alpha
114  ,const DMatrixSlice& gms_rhs1, BLAS_Cpp::Transp trans_rhs1
115  ,BLAS_Cpp::Transp trans_rhs2 ) const;
117  virtual void M_StInvMtM(
118  DMatrix* gm_lhs, value_type alpha
119  ,BLAS_Cpp::Transp trans_rhs1
120  ,const MatrixOpSerial& mwo_rhs2, BLAS_Cpp::Transp trans_rhs2 ) const;
122  virtual void M_StInvMtM(
123  DMatrixSlice* gms_lhs, value_type alpha
124  ,BLAS_Cpp::Transp trans_rhs1
125  ,const MatrixOpSerial& mwo_rhs2, BLAS_Cpp::Transp trans_rhs2 ) const;
127  virtual void M_StMtInvM(
128  DMatrix* gm_lhs, value_type alpha
129  ,const MatrixOpSerial& mwo_rhs1, BLAS_Cpp::Transp trans_rhs1
130  ,BLAS_Cpp::Transp trans_rhs2 ) const;
132  virtual void M_StMtInvM(
133  DMatrixSlice* gms_lhs, value_type alpha
134  ,const MatrixOpSerial& mwo_rhs1, BLAS_Cpp::Transp trans_rhs1
135  ,BLAS_Cpp::Transp trans_rhs2 ) const;
136 
137  // end Level-3 BLAS
139 
142 
144  void V_InvMtV(
145  VectorMutable* v_lhs, BLAS_Cpp::Transp trans_rhs1
146  ,const Vector& v_rhs2) const;
148  void V_InvMtV(
149  VectorMutable* v_lhs, BLAS_Cpp::Transp trans_rhs1
150  ,const SpVectorSlice& sv_rhs2) const;
152  value_type transVtInvMtV(
153  const Vector& v_rhs1
154  ,BLAS_Cpp::Transp trans_rhs2
155  ,const Vector& v_rhs3) const;
157  void M_StInvMtM(
158  MatrixOp* m_lhs, value_type alpha
159  ,BLAS_Cpp::Transp trans_rhs1
160  ,const MatrixOp& mwo_rhs2, BLAS_Cpp::Transp trans_rhs2
161  ) const;
163  void M_StMtInvM(
164  MatrixOp* m_lhs, value_type alpha
165  ,const MatrixOp& mwo_rhs1, BLAS_Cpp::Transp trans_rhs1
166  ,BLAS_Cpp::Transp trans_rhs2
167  ) const;
168 
170 
171 }; // end class MatrixNonsingSerial
172 
179 
182 
184 inline void V_InvMtV(DVector* v_lhs, const MatrixNonsingSerial& M_rhs1
185  , BLAS_Cpp::Transp trans_rhs1, const DVectorSlice& vs_rhs2)
186 {
187  M_rhs1.V_InvMtV(v_lhs,trans_rhs1,vs_rhs2);
188 }
189 
191 inline void V_InvMtV(DVectorSlice* vs_lhs, const MatrixNonsingSerial& M_rhs1
192  , BLAS_Cpp::Transp trans_rhs1, const DVectorSlice& vs_rhs2)
193 {
194  M_rhs1.V_InvMtV(vs_lhs,trans_rhs1,vs_rhs2);
195 }
196 
198 inline void V_InvMtV(DVector* v_lhs, const MatrixNonsingSerial& M_rhs1
199  , BLAS_Cpp::Transp trans_rhs1, const SpVectorSlice& sv_rhs2)
200 {
201  M_rhs1.V_InvMtV(v_lhs,trans_rhs1,sv_rhs2);
202 }
203 
205 inline void V_InvMtV(DVectorSlice* vs_lhs, const MatrixNonsingSerial& M_rhs1
206  , BLAS_Cpp::Transp trans_rhs1, const SpVectorSlice& sv_rhs2)
207 {
208  M_rhs1.V_InvMtV(vs_lhs,trans_rhs1,sv_rhs2);
209 }
210 
212 inline value_type transVtInvMtV(const DVectorSlice& vs_rhs1, const MatrixNonsingSerial& M_rhs2
213  , BLAS_Cpp::Transp trans_rhs2, const DVectorSlice& sv_rhs3)
214 {
215  return M_rhs2.transVtInvMtV(vs_rhs1,trans_rhs2,sv_rhs3);
216 }
217 
219 inline value_type transVtInvMtV(const SpVectorSlice& sv_rhs1, const MatrixNonsingSerial& M_rhs2
220  , BLAS_Cpp::Transp trans_rhs2, const SpVectorSlice& sv_rhs3)
221 {
222  return M_rhs2.transVtInvMtV(sv_rhs1,trans_rhs2,sv_rhs3);
223 }
224 
225 // end Level-2 BLAS
227 
230 
232 inline void M_StInvMtM(
233  DMatrix* gm_lhs, value_type alpha
234  ,const MatrixNonsingSerial& M_rhs1, BLAS_Cpp::Transp trans_rhs1
235  ,const DMatrixSlice& gms_rhs2, BLAS_Cpp::Transp trans_rhs2
236  )
237 {
238  M_rhs1.M_StInvMtM(gm_lhs,alpha,trans_rhs1,gms_rhs2,trans_rhs2);
239 }
240 
242 inline void M_StInvMtM(
243  DMatrixSlice* gms_lhs, value_type alpha
244  ,const MatrixNonsingSerial& M_rhs1, BLAS_Cpp::Transp trans_rhs1
245  ,const DMatrixSlice& gms_rhs2, BLAS_Cpp::Transp trans_rhs2
246  )
247 {
248  M_rhs1.M_StInvMtM(gms_lhs,alpha,trans_rhs1,gms_rhs2,trans_rhs2);
249 }
250 
252 inline void M_StMtInvM(
253  DMatrix* gm_lhs, value_type alpha
254  ,const DMatrixSlice& gms_rhs1, BLAS_Cpp::Transp trans_rhs1
255  ,const MatrixNonsingSerial& M_rhs2, BLAS_Cpp::Transp trans_rhs2
256  )
257 {
258  M_rhs2.M_StMtInvM(gm_lhs,alpha,gms_rhs1,trans_rhs1,trans_rhs2);
259 }
260 
262 inline void M_StMtInvM(
263  DMatrixSlice* gms_lhs, value_type alpha
264  ,const DMatrixSlice& gms_rhs1, BLAS_Cpp::Transp trans_rhs1
265  ,const MatrixNonsingSerial& M_rhs2, BLAS_Cpp::Transp trans_rhs2
266  )
267 {
268  M_rhs2.M_StMtInvM(gms_lhs,alpha,gms_rhs1,trans_rhs1,trans_rhs2);
269 }
270 
272 inline void M_StInvMtM(
273  DMatrix* gm_lhs, value_type alpha
274  ,const MatrixNonsingSerial& M_rhs1, BLAS_Cpp::Transp trans_rhs1
275  ,const MatrixOpSerial& mwo_rhs2, BLAS_Cpp::Transp trans_rhs2
276  )
277 {
278  M_rhs1.M_StInvMtM(gm_lhs,alpha,trans_rhs1,mwo_rhs2,trans_rhs2);
279 }
280 
282 inline void M_StInvMtM(
283  DMatrixSlice* gms_lhs, value_type alpha
284  ,const MatrixNonsingSerial& M_rhs1, BLAS_Cpp::Transp trans_rhs1
285  ,const MatrixOpSerial& mwo_rhs2, BLAS_Cpp::Transp trans_rhs2
286  )
287 {
288  M_rhs1.M_StInvMtM(gms_lhs,alpha,trans_rhs1,mwo_rhs2,trans_rhs2);
289 }
290 
292 inline void M_StMtInvM(
293  DMatrix* gm_lhs, value_type alpha
294  ,const MatrixOpSerial& mwo_rhs1, BLAS_Cpp::Transp trans_rhs1
295  ,const MatrixNonsingSerial& M_rhs2, BLAS_Cpp::Transp trans_rhs2
296  )
297 {
298  M_rhs2.M_StMtInvM(gm_lhs,alpha,mwo_rhs1,trans_rhs1,trans_rhs2);
299 }
300 
302 inline void M_StMtInvM(
303  DMatrixSlice* gms_lhs, value_type alpha
304  ,const MatrixOpSerial& mwo_rhs1, BLAS_Cpp::Transp trans_rhs1
305  ,const MatrixNonsingSerial& M_rhs2, BLAS_Cpp::Transp trans_rhs2
306  )
307 {
308  M_rhs2.M_StMtInvM(gms_lhs,alpha,mwo_rhs1,trans_rhs1,trans_rhs2);
309 }
310 
311 // end Level-3 BLAS
313 
314 // end Inline non-member operation functions
316 
317 } // end namespace AbstractLinAlgPack
318 
319 #endif // SLAP_MATRIX_NONSINGULAR_SERIAL_H
Abstract interface for immutable, finite dimensional, coordinate vectors {abstract}.
void M_StMtInvM(DMatrixSlice *gms_lhs, value_type alpha, const MatrixOpSerial &mwo_rhs1, BLAS_Cpp::Transp trans_rhs1, const MatrixNonsingSerial &M_rhs2, BLAS_Cpp::Transp trans_rhs2)
gms_lhs = alpha * op(mwo_rhs1) * inv(op(M_rhs2)) (left)
Base class for all matrices implemented in a shared memory address space.
virtual void M_StInvMtM(DMatrix *gm_lhs, value_type alpha, BLAS_Cpp::Transp trans_rhs1, const DMatrixSlice &gms_rhs2, BLAS_Cpp::Transp trans_rhs2) const
gm_lhs = alpha * inv(op(M_rhs1)) * op(gms_rhs2) (right)
void V_InvMtV(DVectorSlice *vs_lhs, const MatrixNonsingSerial &M_rhs1, BLAS_Cpp::Transp trans_rhs1, const SpVectorSlice &sv_rhs2)
vs_lhs = inv(op(M_rhs1)) * sv_rhs2
Abstract base class for all AbstractLinAlgPack::MatrixNonsing objects implemented in shared memory sp...
virtual void V_InvMtV(DVector *v_lhs, BLAS_Cpp::Transp trans_rhs1, const DVectorSlice &vs_rhs2) const
v_lhs = inv(op(M_rhs1)) * vs_rhs2
Base class for all matrices that support basic matrix operations.
value_type transVtInvMtV(const SpVectorSlice &sv_rhs1, const MatrixNonsingSerial &M_rhs2, BLAS_Cpp::Transp trans_rhs2, const SpVectorSlice &sv_rhs3)
result = sv_rhs1' * inv(op(M_rhs2)) * sv_rhs3
virtual value_type transVtInvMtV(const DVectorSlice &vs_rhs1, BLAS_Cpp::Transp trans_rhs2, const DVectorSlice &vs_rhs3) const
result = vs_rhs1' * inv(op(M_rhs2)) * vs_rhs3
Abstract interface for mutable coordinate vectors {abstract}.
virtual void M_StMtInvM(DMatrix *gm_lhs, value_type alpha, const DMatrixSlice &gms_rhs1, BLAS_Cpp::Transp trans_rhs1, BLAS_Cpp::Transp trans_rhs2) const
gm_lhs = alpha * op(gms_rhs1) * inv(op(M_rhs2)) (left)
void M_StInvMtM(DMatrixSlice *gms_lhs, value_type alpha, const MatrixNonsingSerial &M_rhs1, BLAS_Cpp::Transp trans_rhs1, const MatrixOpSerial &mwo_rhs2, BLAS_Cpp::Transp trans_rhs2)
gms_lhs = alpha * inv(op(M_rhs1)) * op(mwo_rhs2) (right)
Abstract base class for all nonsingular polymorphic matrices that can solve for linear system with bu...
Transp