|
17 | 17 | */ |
18 | 18 |
|
19 | 19 | /// \file |
| 20 | +/// \brief Functions for caching and calculating Wigner d and D functions |
| 21 | +/// \author Daniel Greenwald, Johannes Rauch |
| 22 | +/// |
| 23 | +/// Calculation is based on M. E. Rose's _Elementary Theory of Angular Momentum_ (1957): |
| 24 | +/// |
| 25 | +/// \f$ D^{J}_{MN}(\alpha, \beta, \gamma) \equiv \exp(-iM\alpha) d^{J}_{MN}(\beta) \exp(-iN\gamma)\f$ |
| 26 | +/// |
| 27 | +/// \f$ d^{J}_{MN}(\beta) \equiv \sum_{K} A (-)^{K+M-N} B_{K} (\cos\frac{\beta}{2})^{2J+N-M-2K} (\sin\frac{\beta}{2})^{M-N+2K} \f$ |
| 28 | +/// |
| 29 | +/// with \f$ A \equiv (J+N)!(J-N)!(J+M)!(J-M)!\f$ |
| 30 | +/// and \f$ B_{K} \equiv (J-M-K)!((J+N-K)!(K+M-N)!K!\f$. |
| 31 | +/// |
| 32 | +/// The limits on K are governed by the factorials, |
| 33 | +/// since \f$ (1 / X!) = 0 \f$ if \f$ X < 0\f$. |
| 34 | +/// |
| 35 | +/// We only cache matrix elements with M in [-J, J] and N in [-J, min(0, m)]. |
| 36 | +/// This amounts to the lower triangle and the diagonal without the bottom right corner. |
| 37 | +/// The uncached matrix elements are given by the by the symmetries |
| 38 | +/// - \f$ d^{J}_{MN}(\beta) = (-)^(M-N) d^{J}_{NM}(\beta)\f$ |
| 39 | +/// - \f$ d^{J}_{MN}(\beta) = (-)^{M-N) d^{J}_{-N-M}(\beta)\f$ |
20 | 40 |
|
21 | 41 | #ifndef yap_WignerD_h |
22 | 42 | #define yap_WignerD_h |
23 | 43 |
|
| 44 | +#include "Constants.h" |
| 45 | + |
24 | 46 | #include <complex> |
25 | | -#include <vector> |
26 | 47 |
|
27 | 48 | namespace yap { |
28 | 49 |
|
29 | | -/// Wigner d-function d^j_{two_m two_n}(theta) |
30 | | -double dFunction(const int two_j, const int two_m, const int two_n, const double theta); |
31 | | - |
32 | | -/// spherical harmonics Y_l^{two_m}(theta, phi) |
33 | | -std::complex<double> sphericalHarmonic(const int two_l, const int two_m, const double theta, const double phi); |
34 | | - |
35 | | -/// Wigner D-function D^j_{two_m two_n}(alpha, beta, gamma) |
36 | | -std::complex<double> DFunction(const int two_j, const int two_m, const int two_n, const double alpha, const double beta, const double gamma); |
37 | | - |
38 | | -/// complex conjugate of Wigner D-function D^j_{two_m two_n}(alpha, beta, gamma) |
39 | | -std::complex<double> DFunctionConj(const int two_j, const int two_m, const int two_n, const double alpha, const double beta, const double gamma); |
40 | | - |
41 | | - |
42 | | -/// \class dFunctionCached |
43 | | -/// \brief cached Wigner d and D functions |
44 | | -/// \author Daniel Greenwald, Johannes Rauch |
45 | | -/// This class has been copied from rootpwa and modified |
46 | | -class dFunctionCached |
47 | | -{ |
48 | | - |
49 | | -public: |
50 | | - |
51 | | - struct cacheEntryType { |
52 | | - cacheEntryType() |
53 | | - : constTerm(0) |
54 | | - { } |
| 50 | +/// Wigner D-function \f$ D^{J}_{M N}(\alpha, \beta, \gamma) \f$ |
| 51 | +std::complex<double> DFunction(unsigned char twoJ, char twoM, char twoN, double alpha, double beta, double gamma); |
55 | 52 |
|
56 | | - double constTerm; |
57 | | - std::vector<int> kmn1; |
58 | | - std::vector<int> jmnk; |
59 | | - std::vector<double> factor; |
60 | | - }; |
| 53 | +/// \return Wigner d-function \f$ d^{J}_{M N}(\beta) \f$ |
| 54 | +/// \param twoJ twice the total spin of system |
| 55 | +/// \param twoM twice the first spin projection |
| 56 | +/// \param twoN twice the second spin projection |
| 57 | +/// \param beta rotation angle |
| 58 | +double dFunction(unsigned char twoJ, char twoM, char twoN, double beta); |
61 | 59 |
|
62 | | - typedef std::vector<std::vector<std::vector<cacheEntryType> > > cacheType; |
63 | 60 |
|
| 61 | +namespace dMatrix { |
64 | 62 |
|
65 | | - /// get singleton instance |
66 | | - static dFunctionCached& instance() |
67 | | - { return _instance; } |
| 63 | +/// Cache d-matrix elements for representation of spin J |
| 64 | +/// \param twoJ twice the spin of the representation |
| 65 | +void cache(unsigned char twoJ); |
68 | 66 |
|
69 | | - /// \return d^j_{two_m two_n}(theta) |
70 | | - double operator ()(const int two_j, const int two_m, const int two_n, const double theta); |
71 | | - |
72 | | - /// \return caching flag |
73 | | - static bool useCache() |
74 | | - { return _useCache; } |
75 | | - |
76 | | - /// sets caching flag |
77 | | - static void setUseCache(const bool useCache = true) |
78 | | - { _useCache = useCache; } |
79 | | - |
80 | | - /// \return cache size in bytes |
81 | | - static unsigned int cacheSize(); |
82 | | - |
83 | | - |
84 | | -private: |
85 | | - |
86 | | - dFunctionCached(); |
87 | | - ~dFunctionCached(); |
88 | | - dFunctionCached(const dFunctionCached&); |
89 | | - dFunctionCached& operator =(const dFunctionCached&); |
90 | | - |
91 | | - static dFunctionCached _instance; ///< singleton instance |
92 | | - static bool _useCache; ///< if set to true, cache is used |
93 | | - |
94 | | - static const unsigned int _maxJ = 21; ///< maximum allowed angular momentum * 2 + 1 |
95 | | - static cacheEntryType* _cache[_maxJ][_maxJ + 1][_maxJ + 1]; ///< cache for intermediate terms [two_j][two_m][two_n] |
96 | | -}; |
| 67 | +/// \return cache size in bytes |
| 68 | +unsigned int cacheSize(); |
97 | 69 |
|
| 70 | +} |
98 | 71 |
|
99 | 72 | } |
100 | 73 |
|
|
0 commit comments