5#ifndef GEMMI_CELLRED_HPP_
6#define GEMMI_CELLRED_HPP_
31 :
A(
m.column_dot(0,0)),
34 xi(2 *
m.column_dot(1,2)),
35 eta(2 *
m.column_dot(0,2)),
36 zeta(2 *
m.column_dot(0,1)) {}
55 double a = std::sqrt(
A);
56 double b = std::sqrt(
B);
57 double c = std::sqrt(
C);
59 deg(std::acos(
xi/(2*b*c))),
60 deg(std::acos(
eta/(2*a*c))),
69 return A <=
B &&
B <=
C &&
70 (
A !=
B || std::abs(
xi) <= std::abs(
eta)) &&
71 (
B !=
C || std::abs(
eta) <= std::abs(
zeta)) &&
91 swap_columns_and_negate(0, 1);
99 swap_columns_and_negate(1, 2);
110 if (
sgn *
xi < -
eps) negate_column(0);
114 negate_column(std::fabs(
zeta) <=
eps ? 2 :
115 std::fabs(
eta) <=
eps ? 1 : 0);
125 if (std::abs(
xi) >
B) {
126 double j = std::floor(0.5*
xi/
B + 0.5);
130 }
else if (std::abs(
eta) >
A) {
131 double j = std::floor(0.5*
eta/
A + 0.5);
135 }
else if (std::abs(
zeta) >
A) {
136 double j = std::floor(0.5*
zeta/
A + 0.5);
141 double j = std::floor(0.5 * (
xi +
eta) / (
A +
B +
zeta) + 0.5);
162 double sum = std::sqrt(
A) + std::sqrt(
B) + std::sqrt(
C);
163 if (std::abs(sum -
prev_sum) < sum * 1
e-6) {
189 add_column(1, 2, -
int(
sign_xi));
239 void swap_columns_and_negate(
int i,
int j) {
246 void negate_column(
int i) {
250 void add_column(
int pos,
int dest,
int sign) {
252 r[dest] += sign * r[pos];
267 std::array<double,6>
s;
273 for (
int i = 0;
i < 3; ++
i)
275 b[3]= -b[0] - b[1] - b[2];
276 s[0] = b[1].dot(b[2]);
277 s[1] = b[0].dot(b[2]);
278 s[2] = b[0].dot(b[1]);
279 s[3] = b[0].dot(b[3]);
280 s[4] = b[1].dot(b[3]);
281 s[5] = b[2].dot(b[3]);
291 return -2 * (
s[0] +
s[1] +
s[2] +
s[3] +
s[4] +
s[5]);
295 return std::all_of(
s.begin(),
s.end(), [
eps](
double x) { return x <= eps; });
301 const int table[6][5] = {
314 for (
int i = 0;
i < 6; ++
i)
343 return {-
s[1]-
s[2]-
s[3], -
s[0]-
s[2]-
s[4], -
s[0]-
s[1]-
s[5], 2*
s[0], 2*
s[1], 2*
s[2]};
352 s[1]+
s[2]+
s[3],
s[0]+
s[2]+
s[4],
s[0]+
s[1]+
s[5],
s[3]+
s[4]+
s[5]
356 for (
int i = 0;
i < 3; ++
i)
361 std::swap(
s[1],
s[5]);
362 std::swap(
s[2],
s[4]);
365 std::swap(
s[0],
s[5]);
366 std::swap(
s[2],
s[3]);
369 std::swap(
s[0],
s[4]);
370 std::swap(
s[1],
s[3]);
376 std::swap(
s[0],
s[1]);
377 std::swap(
s[3],
s[4]);
381 std::swap(
s[1],
s[2]);
382 std::swap(
s[4],
s[5]);
386 std::swap(
s[0],
s[1]);
387 std::swap(
s[3],
s[4]);
399 double s0 = 0.5 *
xi;
400 double s1 = 0.5 *
eta;
constexpr double deg(double angle)
Op::Rot centred_to_primitive(char centring_type)
GruberVector(const std::array< double, 6 > &g6)
void set_change_of_basis(const Op &op)
void normalize(double eps=1e-9)
int niggli_reduce(double epsilon=1e-9, int iteration_limit=100)
bool niggli_step(double epsilon=1e-9)
GruberVector(const UnitCell &u, const SpaceGroup *sg, bool track_change_of_basis=false)
SellingVector selling() const
bool is_buerger(double epsilon=1e-9) const
std::unique_ptr< Op > change_of_basis
bool is_normalized() const
bool is_niggli(double epsilon=1e-9) const
GruberVector(const UnitCell &u, char centring, bool track_change_of_basis=false)
UnitCell get_cell() const
GruberVector(const Mat33 &m)
std::array< double, 6 > parameters() const
std::array< double, 6 > cell_parameters() const
Vec3 column_copy(int i) const
std::array< double, 6 > s
bool is_reduced(double eps=1e-9) const
bool reduce_step(double eps=1e-9)
std::array< double, 6 > g6_parameters() const
SellingVector(const Mat33 &orth)
SellingVector(const UnitCell &u, char centring)
UnitCell get_cell() const
std::array< double, 6 > cell_parameters() const
int reduce(double eps=1e-9, int iteration_limit=100)
SellingVector(const UnitCell &u, const SpaceGroup *sg)
void sort(double eps=1e-9)
SellingVector(const std::array< double, 6 > &s_)
double sum_b_squared() const
GruberVector gruber() const