sr1.hpp 1.17 KB
Newer Older
1
// Software License for MTL
2
3
//
// Copyright (c) 2007 The Trustees of Indiana University.
4
//               2008 Dresden University of Technology and the Trustees of Indiana University.
5
//               2010 SimuNova UG (haftungsbeschränkt), www.simunova.com.
6
7
// All rights reserved.
// Authors: Peter Gottschling and Andrew Lumsdaine
8
//
9
// This file is part of the Matrix Template Library
10
//
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
// See also license.mtl.txt in the distribution.

#ifndef ITL_SR1_INCLUDE
#define ITL_SR1_INCLUDE

#include <boost/numeric/mtl/concept/collection.hpp>
#include <boost/numeric/mtl/utility/exception.hpp>
#include <boost/numeric/itl/utility/exception.hpp>

namespace itl {

/// Update of Hessian matrix for e.g. Quasi-Newton by SR1 formulas
struct sr1
{
    /// \f$ H_{k+1}=B_{k+1}^{-1}=H_k+\frac{(s_k-H_k\cdot y_k)(s_k-H_k\cdot y_k)^T}{(s_k-H_k\cdot y_k)^T\cdot y_k} \f$
    template <typename Matrix, typename Vector>
    void operator() (Matrix& H, const Vector& y, const Vector& s)
    {
	assert(num_rows(H) == num_cols(H));
	Vector     d(s - H * y);
	MTL_THROW_IF(dot(d, y) == 0.0, unexpected_orthogonality());
	H+= 1 / dot(d, y) * d * trans(d);
   }
34
};
35
36
37
38
39
40



} // namespace itl

#endif // ITL_SR1_INCLUDE