Householder.h
Go to the documentation of this file.
1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2010 Benoit Jacob <jacob.benoit.1@gmail.com>
5 // Copyright (C) 2009 Gael Guennebaud <gael.guennebaud@inria.fr>
6 //
7 // Eigen is free software; you can redistribute it and/or
8 // modify it under the terms of the GNU Lesser General Public
9 // License as published by the Free Software Foundation; either
10 // version 3 of the License, or (at your option) any later version.
11 //
12 // Alternatively, you can redistribute it and/or
13 // modify it under the terms of the GNU General Public License as
14 // published by the Free Software Foundation; either version 2 of
15 // the License, or (at your option) any later version.
16 //
17 // Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
18 // WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
19 // FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the
20 // GNU General Public License for more details.
21 //
22 // You should have received a copy of the GNU Lesser General Public
23 // License and a copy of the GNU General Public License along with
24 // Eigen. If not, see <http://www.gnu.org/licenses/>.
25 
26 #ifndef EIGEN_HOUSEHOLDER_H
27 #define EIGEN_HOUSEHOLDER_H
28 
29 namespace Eigen {
30 
31 namespace internal {
32 template<int n> struct decrement_size
33 {
34  enum {
35  ret = n==Dynamic ? n : n-1
36  };
37 };
38 }
39 
56 template<typename Derived>
58 {
60  makeHouseholder(essentialPart, tau, beta);
61 }
62 
78 template<typename Derived>
79 template<typename EssentialPart>
81  EssentialPart& essential,
82  Scalar& tau,
83  RealScalar& beta) const
84 {
85  EIGEN_STATIC_ASSERT_VECTOR_ONLY(EssentialPart)
87 
88  RealScalar tailSqNorm = size()==1 ? RealScalar(0) : tail.squaredNorm();
89  Scalar c0 = coeff(0);
90 
91  if(tailSqNorm == RealScalar(0) && internal::imag(c0)==RealScalar(0))
92  {
93  tau = RealScalar(0);
94  beta = internal::real(c0);
95  essential.setZero();
96  }
97  else
98  {
99  beta = internal::sqrt(internal::abs2(c0) + tailSqNorm);
100  if (internal::real(c0)>=RealScalar(0))
101  beta = -beta;
102  essential = tail / (c0 - beta);
103  tau = internal::conj((beta - c0) / beta);
104  }
105 }
106 
122 template<typename Derived>
123 template<typename EssentialPart>
125  const EssentialPart& essential,
126  const Scalar& tau,
127  Scalar* workspace)
128 {
129  if(rows() == 1)
130  {
131  *this *= Scalar(1)-tau;
132  }
133  else
134  {
137  tmp.noalias() = essential.adjoint() * bottom;
138  tmp += this->row(0);
139  this->row(0) -= tau * tmp;
140  bottom.noalias() -= tau * essential * tmp;
141  }
142 }
143 
159 template<typename Derived>
160 template<typename EssentialPart>
162  const EssentialPart& essential,
163  const Scalar& tau,
164  Scalar* workspace)
165 {
166  if(cols() == 1)
167  {
168  *this *= Scalar(1)-tau;
169  }
170  else
171  {
174  tmp.noalias() = right * essential.conjugate();
175  tmp += this->col(0);
176  this->col(0) -= tau * tmp;
177  right.noalias() -= tau * tmp * essential.transpose();
178  }
179 }
180 
181 } // end namespace Eigen
182 
183 #endif // EIGEN_HOUSEHOLDER_H