15 #include <Eigen/Dense>
21 static constexpr
int LTmn = LogArea::Tmn::id;
24 using Eigen::Matrix4d;
25 using Eigen::Vector4d;
49 logg[
LTmn].debug(
"Looking for Landau frame for T_{mu}^{nu} ", A);
50 Eigen::EigenSolver<Matrix4d> es(A);
60 Vector4d eig_im = es.eigenvalues().imag();
61 Vector4d eig_re = es.eigenvalues().real();
62 for (
size_t i = 0; i < 4; i++) {
71 assert(eig_re(0) >= eig_re(1));
72 assert(eig_re(0) >= eig_re(2));
73 assert(eig_re(0) >= eig_re(3));
74 logg[
LTmn].debug(
"eigenvalues: ", eig_re);
77 Vector4d tmp = es.eigenvectors().col(0).real();
84 const double u_sqr = u.
sqr();
86 u /= std::sqrt(u_sqr);
89 "Landau frame is not defined.",
" Eigen vector", u,
" of ", A,
90 " is not time-like and",
91 " cannot be 4-velocity. This may happen if energy-momentum",
92 " tensor was constructed for a massless particle.");
99 using Eigen::Matrix4d;
111 L << u[0], u[1], u[2], u[3],
112 u[1], u[1] * tmp.
x1() + 1.0, u[2] * tmp.
x1(), u[3] * tmp.
x1(),
113 u[2], u[1] * tmp.
x2(), u[2] * tmp.
x2() + 1.0, u[3] * tmp.
x2(),
114 u[3], u[1] * tmp.
x3(), u[2] * tmp.
x3(), u[3] * tmp.
x3() + 1.0;
120 R(1, 1),
R(1, 2),
R(1, 3),
132 Tmn_[4] += mom[1] * tmp.
x1();
133 Tmn_[5] += mom[1] * tmp.
x2();
134 Tmn_[6] += mom[1] * tmp.
x3();
135 Tmn_[7] += mom[2] * tmp.
x2();
136 Tmn_[8] += mom[2] * tmp.
x3();
137 Tmn_[9] += mom[3] * tmp.
x3();
146 for (
size_t mu = 0; mu < 4; mu++) {
147 for (
size_t nu = 0; nu < 4; nu++) {
148 out << std::setprecision(3) << std::setw(12) << std::fixed