details.hpp 1.74 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
/******************************************************************************
 *
 * AMDiS - Adaptive multidimensional simulations
 *
 * Copyright (C) 2013 Dresden University of Technology. All Rights Reserved.
 * Web: https://fusionforge.zih.tu-dresden.de/projects/amdis
 *
 * Authors: 
 * Simon Vey, Thomas Witkowski, Andreas Naumann, Simon Praetorius, et al.
 *
 * This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
 * WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
 *
 *
 * This file is part of AMDiS
 *
 * See also license.opensource.txt in the distribution.
 * 
 ******************************************************************************/

// Written by Simon Praetorius


#ifndef ITL_DETAIL_INCLUDE
#define ITL_DETAIL_INCLUDE

#include <boost/math/special_functions/sign.hpp>

namespace itl {
  
  namespace details {
    
    /// Compute the Givens rotation matrix parameters for a and b.
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
//     template<typename T>
//     void rotmat(const T& a, const T& b , T& c, T& s)
//     {
//       using std::abs; using std::sqrt; using mtl::conj;
//       
//       const T zero = math::zero(T());
//       if (a == zero) {
// 	c = 0.0;
// 	s = 1.0;
//       } else {
// 	double temp = abs(a) / sqrt( conj(a)*a + conj(b)*b );
// 	c = temp;
// 	s = temp * (b / a);
//       }
//     }
    
    inline void rotmat(const double& a, const double& b , double& c, double& s)
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
    {
      using std::abs; using std::sqrt;
      if ( b == 0.0 ) {
	  c = 1.0;
	  s = 0.0;
      } else if ( abs(b) > abs(a) ) {
	  double temp = a / b;
	  s = 1.0 / sqrt( 1.0 + temp*temp );
	  c = temp * s;
      } else {
	  double temp = b / a;
	  c = 1.0 / sqrt( 1.0 + temp*temp );
	  s = temp * c;
      }
    }
  }
}

#endif // ITL_DETAIL_INCLUDE