1 ///////////////////////////////////////////////////////////////////////////////
2 //                                                                           //
3 // The Template Matrix/Vector Library for C++ was created by Mike Jarvis     //
4 // Copyright (C) 1998 - 2016                                                 //
5 // All rights reserved                                                       //
6 //                                                                           //
7 // The project is hosted at https://code.google.com/p/tmv-cpp/               //
8 // where you can find the current version and current documention.           //
9 //                                                                           //
10 // For concerns or problems with the software, Mike may be contacted at      //
11 // mike_jarvis17 [at] gmail.                                                 //
12 //                                                                           //
13 // This software is licensed under a FreeBSD license.  The file              //
14 // TMV_LICENSE should have bee included with this distribution.              //
15 // It not, you can get a copy from https://code.google.com/p/tmv-cpp/.       //
16 //                                                                           //
17 // Essentially, you can use this software however you want provided that     //
18 // you include the TMV_LICENSE file in any distribution that uses it.        //
19 //                                                                           //
20 ///////////////////////////////////////////////////////////////////////////////
21 
22 #ifndef TMV_BandLUDiv_H
23 #define TMV_BandLUDiv_H
24 
25 #include "tmv/TMV_BaseBandMatrix.h"
26 #include "tmv/TMV_BaseTriMatrix.h"
27 
28 namespace tmv {
29 
30 
31     template <typename T, typename T1>
32     void TriLDivEq(
33         const GenBandMatrix<T1>& A, MatrixView<T> B, DiagType dt);
34     template <typename T, typename T1>
35     void TriLDivEq(
36         const GenBandMatrix<T1>& A, VectorView<T> v, DiagType dt);
37     // Solve A X = B  where A is upper or lower band triangular
38 
39     template <typename T, typename T1>
40     void LU_PackedPL_Unpack(
41         const GenBandMatrix<T1>& LUx, const ptrdiff_t* p,
42         LowerTriMatrixView<T> m);
43 
44     template <typename T, typename T1>
45     void LU_PackedPL_LDivEq(
46         const GenBandMatrix<T1>& LUx, const ptrdiff_t* p, MatrixView<T> m);
47     template <typename T, typename T1>
48     void LU_PackedPL_RDivEq(
49         const GenBandMatrix<T1>& LUx, const ptrdiff_t* p, MatrixView<T> m);
50 
51     template <typename T, typename T1>
52     void LU_LDivEq(
53         const GenBandMatrix<T1>& LUx, const ptrdiff_t* p, MatrixView<T> m);
54     template <typename T, typename T1>
55     void LU_RDivEq(
56         const GenBandMatrix<T1>& LUx, const ptrdiff_t* p, MatrixView<T> m);
57 
58     template <typename T, typename T1>
59     void LU_Inverse(
60         const GenBandMatrix<T1>& LUx, const ptrdiff_t* p, MatrixView<T> m);
61 
62     template <typename T>
63     void TriInverse(UpperTriMatrixView<T> U, ptrdiff_t nhi);
64 
65 #ifndef NOTHROW
66     template <typename T>
67     class SingularBandLU : public Singular
68     {
69     public:
70 
71         BandMatrix<T> A;
72 
SingularBandLU(const GenBandMatrix<T> & _A)73         SingularBandLU(const GenBandMatrix<T>& _A) :
74             Singular("BandMatrix."), A(_A) {}
throw()75         ~SingularBandLU() throw() {}
76 
write(std::ostream & os)77         void write(std::ostream& os) const throw()
78         {
79             Singular::write(os);
80             os<<"In LU Decomposed form, the matrix is \n"<<A<<std::endl;
81             os<<"ie. U = "<<A.upperBand()<<std::endl;
82         }
83     };
84 
85     template <typename T>
86     class SingularBandLU2 : public SingularBandLU<T>
87     {
88     public:
89 
90         BandMatrix<T> A0;
91 
SingularBandLU2(const GenBandMatrix<T> & _A,const GenBandMatrix<T> & _A0)92         SingularBandLU2(
93             const GenBandMatrix<T>& _A, const GenBandMatrix<T>& _A0) :
94             SingularBandLU<T>(_A), A0(_A0) {}
throw()95         ~SingularBandLU2() throw() {}
96 
write(std::ostream & os)97         void write(std::ostream& os) const throw()
98         {
99             SingularBandLU<T>::write(os);
100             os<<"The original BandMatrix is \n"<<A0<<std::endl;
101         }
102     };
103 #endif
104 
105 
106     // Specialize disallowed complex combinations:
107 #define CT std::complex<T>
108     template <typename T>
TriLDivEq(const GenBandMatrix<CT> &,MatrixView<T>,DiagType)109     inline void TriLDivEq(
110         const GenBandMatrix<CT>& , MatrixView<T> , DiagType )
111     { TMVAssert(TMV_FALSE); }
112     template <typename T>
TriLDivEq(const GenBandMatrix<CT> &,VectorView<T>,DiagType)113     inline void TriLDivEq(
114         const GenBandMatrix<CT>& , VectorView<T> , DiagType )
115     { TMVAssert(TMV_FALSE); }
116 
117     template <typename T>
LU_PackedPL_Unpack(const GenBandMatrix<CT> &,const ptrdiff_t *,LowerTriMatrixView<T>)118     inline void LU_PackedPL_Unpack(
119         const GenBandMatrix<CT>& , const ptrdiff_t* ,
120         LowerTriMatrixView<T> )
121     { TMVAssert(TMV_FALSE); }
122 
123     template <typename T>
LU_PackedPL_LDivEq(const GenBandMatrix<CT> &,const ptrdiff_t *,MatrixView<T>)124     inline void LU_PackedPL_LDivEq(
125         const GenBandMatrix<CT>& , const ptrdiff_t* , MatrixView<T> )
126     { TMVAssert(TMV_FALSE); }
127     template <typename T>
LU_PackedPL_RDivEq(const GenBandMatrix<CT> &,const ptrdiff_t *,MatrixView<T>)128     inline void LU_PackedPL_RDivEq(
129         const GenBandMatrix<CT>& , const ptrdiff_t* , MatrixView<T> )
130     { TMVAssert(TMV_FALSE); }
131 
132     template <typename T>
LU_LDivEq(const GenBandMatrix<CT> &,const ptrdiff_t *,MatrixView<T>)133     inline void LU_LDivEq(
134         const GenBandMatrix<CT>& , const ptrdiff_t* , MatrixView<T> )
135     { TMVAssert(TMV_FALSE); }
136     template <typename T>
LU_RDivEq(const GenBandMatrix<CT> &,const ptrdiff_t *,MatrixView<T>)137     inline void LU_RDivEq(
138         const GenBandMatrix<CT>& , const ptrdiff_t* , MatrixView<T> )
139     { TMVAssert(TMV_FALSE); }
140 
141     template <typename T>
LU_Inverse(const GenBandMatrix<CT> &,const ptrdiff_t *,MatrixView<T>)142     inline void LU_Inverse(
143         const GenBandMatrix<CT>& , const ptrdiff_t* , MatrixView<T> )
144     { TMVAssert(TMV_FALSE); }
145 
146 #undef CT
147 
148 }
149 
150 #endif
151