leonarski_f
527a5187f4
ITA character 43 was absent from the table. It is the type-II reduced form of a centred monoclinic lattice with no length equality - mC, mI, mA and mF are one Bravais lattice in four settings, and this row names the one whose conventional cell comes out I-centred. With the row missing such a cell reaches character 44 and is reported as triclinic, losing its centring outright. It accounted for 30 of the 31 demotions left after the two fixes before this one, and for 4.1% of random centred-monoclinic lattices. Both of its conditions are equalities on scalar products - International Tables gives them as 2|D+E+F| = A+B and |2D+F| = B - so the three angle tests are vacuous for it and the second equality is what selects it. cond_2DF was declared in the character struct and never tested against anything, because until now no row used it. Audited over 14000 exact lattices, the row fires 35 times: 32 are centred monoclinic lattices it recovers correctly and 3 are triclinic cells it promotes, all three of which Le Page promotes at the same tolerance. No over-call is attributable to it. The row is taken from International Tables rather than derived - no single integer matrix of determinant 2 covers more than 68% of these cells - and the test checks the conventional cell it produces has two right angles and twice the primitive volume. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01EFEJG6WBQv8th4UJFNe53N (cherry picked from commit 4c2baf0cc3b2177196a9448df06e24b5446de96a)
577 lines
22 KiB
C++
577 lines
22 KiB
C++
// SPDX-FileCopyrightText: 2025 Filip Leonarski, Paul Scherrer Institute <filip.leonarski@psi.ch>
|
|
// SPDX-License-Identifier: GPL-3.0-only
|
|
|
|
#include <catch2/catch_all.hpp>
|
|
#include "../common/CrystalLattice.h"
|
|
#include "../common/Coord.h"
|
|
#include "../common/UnitCell.h"
|
|
#include "../image_analysis/lattice_search/LatticeSearch.h"
|
|
#include "gemmi/symmetry.hpp"
|
|
#include <cmath>
|
|
|
|
// Helper: check near-equality of unit cell parameters
|
|
static void check_uc(const UnitCell& uc, double a, double b, double c,
|
|
double alpha, double beta, double gamma,
|
|
double eps_len = 1e-6, double eps_ang = 1e-4) {
|
|
CHECK(uc.a == Catch::Approx(a).margin(eps_len));
|
|
CHECK(uc.b == Catch::Approx(b).margin(eps_len));
|
|
CHECK(uc.c == Catch::Approx(c).margin(eps_len));
|
|
CHECK(uc.alpha == Catch::Approx(alpha).margin(eps_ang));
|
|
CHECK(uc.beta == Catch::Approx(beta ).margin(eps_ang));
|
|
CHECK(uc.gamma == Catch::Approx(gamma).margin(eps_ang));
|
|
}
|
|
|
|
|
|
TEST_CASE("LatticeSearch - cubic I") {
|
|
// Build a body-centered cubic cell with a=40:
|
|
// primitive basis vectors (conventional I cubic primitive):
|
|
// p1 = (0, a/2, a/2), p2 = (a/2, 0, a/2), p3 = (a/2, a/2, 0)
|
|
const double a = 40.0;
|
|
CrystalLattice L(
|
|
Coord(a, 0, 0),
|
|
Coord(0, a, 0),
|
|
Coord(0, 0, a)
|
|
);
|
|
L = L.ToPrimitive('I');
|
|
|
|
auto res = LatticeSearch(L, 1e-6);
|
|
|
|
CHECK(res.system == gemmi::CrystalSystem::Cubic);
|
|
CHECK(res.centering == 'I');
|
|
|
|
// Conventional cubic I should have equal edges and 90° angles
|
|
auto uc = res.conventional.GetUnitCell();
|
|
CHECK(uc.a == Catch::Approx( a )); // In this construction, conventional a matches given a
|
|
CHECK(uc.b == Catch::Approx( a ));
|
|
CHECK(uc.c == Catch::Approx( a ));
|
|
CHECK(uc.alpha == Catch::Approx(90.0));
|
|
CHECK(uc.beta == Catch::Approx(90.0));
|
|
CHECK(uc.gamma == Catch::Approx(90.0));
|
|
}
|
|
|
|
TEST_CASE("LatticeSearch - cubic F") {
|
|
// Build a body-centered cubic cell with a=40:
|
|
// primitive basis vectors (conventional I cubic primitive):
|
|
// p1 = (0, a/2, a/2), p2 = (a/2, 0, a/2), p3 = (a/2, a/2, 0)
|
|
const double a = 40.0;
|
|
CrystalLattice L(
|
|
Coord(a, 0, 0),
|
|
Coord(0, a, 0),
|
|
Coord(0, 0, a)
|
|
);
|
|
L = L.ToPrimitive('F');
|
|
|
|
auto res = LatticeSearch(L, 1e-6);
|
|
|
|
CHECK(res.system == gemmi::CrystalSystem::Cubic);
|
|
CHECK(res.centering == 'F');
|
|
|
|
// Conventional cubic I should have equal edges and 90° angles
|
|
auto uc = res.conventional.GetUnitCell();
|
|
CHECK(uc.a == Catch::Approx( a )); // In this construction, conventional a matches given a
|
|
CHECK(uc.b == Catch::Approx( a ));
|
|
CHECK(uc.c == Catch::Approx( a ));
|
|
CHECK(uc.alpha == Catch::Approx(90.0));
|
|
CHECK(uc.beta == Catch::Approx(90.0));
|
|
CHECK(uc.gamma == Catch::Approx(90.0));
|
|
}
|
|
|
|
TEST_CASE("LatticeSearch - cubic P") {
|
|
// Simple cubic P, a=30
|
|
const double a = 30.0;
|
|
CrystalLattice L(
|
|
Coord(a,0,0),
|
|
Coord(0,a,0),
|
|
Coord(0,0,a)
|
|
);
|
|
|
|
auto res = LatticeSearch(L, 1e-6);
|
|
|
|
CHECK(res.system == gemmi::CrystalSystem::Cubic);
|
|
CHECK(res.centering == 'P');
|
|
|
|
auto uc = res.conventional.GetUnitCell();
|
|
check_uc(uc, a, a, a, 90.0, 90.0, 90.0, 1e-6, 1e-4);
|
|
}
|
|
|
|
TEST_CASE("LatticeSearch - tetragonal I") {
|
|
// Build a body-centered cubic cell with a=40:
|
|
// primitive basis vectors (conventional I cubic primitive):
|
|
// p1 = (0, a/2, a/2), p2 = (a/2, 0, a/2), p3 = (a/2, a/2, 0)
|
|
const double a = 40.0;
|
|
const double b = 34.0;
|
|
CrystalLattice L(
|
|
Coord(a, 0, 0),
|
|
Coord(0, a, 0),
|
|
Coord(0, 0, b)
|
|
);
|
|
L = L.ToPrimitive('I');
|
|
|
|
auto res = LatticeSearch(L, 1e-6);
|
|
|
|
CHECK(res.system == gemmi::CrystalSystem::Tetragonal);
|
|
CHECK(res.centering == 'I');
|
|
|
|
// Conventional cubic I should have equal edges and 90° angles
|
|
auto uc = res.conventional.GetUnitCell();
|
|
CHECK(uc.a == Catch::Approx( a )); // In this construction, conventional a matches given a
|
|
CHECK(uc.b == Catch::Approx( a ));
|
|
CHECK(uc.c == Catch::Approx( b ));
|
|
CHECK(uc.alpha == Catch::Approx(90.0));
|
|
CHECK(uc.beta == Catch::Approx(90.0));
|
|
CHECK(uc.gamma == Catch::Approx(90.0));
|
|
}
|
|
|
|
TEST_CASE("LatticeSearch - tetragonal I - v2") {
|
|
// Build a body-centered cubic cell with a=40:
|
|
// primitive basis vectors (conventional I cubic primitive):
|
|
// p1 = (0, a/2, a/2), p2 = (a/2, 0, a/2), p3 = (a/2, a/2, 0)
|
|
const double a = 40.0;
|
|
const double b = 54.0;
|
|
CrystalLattice L(
|
|
Coord(a, 0, 0),
|
|
Coord(0, a, 0),
|
|
Coord(0, 0, b)
|
|
);
|
|
L = L.ToPrimitive('I');
|
|
|
|
auto res = LatticeSearch(L, 1e-6);
|
|
|
|
CHECK(res.system == gemmi::CrystalSystem::Tetragonal);
|
|
CHECK(res.centering == 'I');
|
|
|
|
// Conventional cubic I should have equal edges and 90° angles
|
|
auto uc = res.conventional.GetUnitCell();
|
|
CHECK(uc.a == Catch::Approx( a )); // In this construction, conventional a matches given a
|
|
CHECK(uc.b == Catch::Approx( a ));
|
|
CHECK(uc.c == Catch::Approx( b ));
|
|
CHECK(uc.alpha == Catch::Approx(90.0));
|
|
CHECK(uc.beta == Catch::Approx(90.0));
|
|
CHECK(uc.gamma == Catch::Approx(90.0));
|
|
}
|
|
|
|
// Tetragonal P: a=b!=c, all angles 90, P-centering
|
|
TEST_CASE("LatticeSearch - tetragonal P") {
|
|
const double a = 37.0, c = 59.0;
|
|
CrystalLattice L(
|
|
Coord(a,0,0),
|
|
Coord(0,a,0),
|
|
Coord(0,0,c)
|
|
);
|
|
|
|
auto res = LatticeSearch(L, 1e-6);
|
|
|
|
CHECK(res.system == gemmi::CrystalSystem::Tetragonal);
|
|
CHECK(res.centering == 'P');
|
|
|
|
auto uc = res.conventional.GetUnitCell();
|
|
check_uc(uc, a, a, c, 90.0, 90.0, 90.0, 1e-2, 1e-2);
|
|
}
|
|
|
|
// Orthorhombic F: all angles 90, unequal edges, F-centering
|
|
TEST_CASE("LatticeSearch - orthorhombic F") {
|
|
const double a = 35.0, b = 41.0, c = 57.0;
|
|
CrystalLattice conv(a,b,c, 90.0,90.0,90.0);
|
|
|
|
CrystalLattice L = conv.ToPrimitive('F');
|
|
|
|
auto res = LatticeSearch(L, 1e-6);
|
|
|
|
CHECK(res.system == gemmi::CrystalSystem::Orthorhombic);
|
|
CHECK(res.centering == 'F');
|
|
|
|
auto uc = res.conventional.GetUnitCell();
|
|
check_uc(uc, a, b, c, 90.0, 90.0, 90.0, 1e-1, 1e-2);
|
|
}
|
|
|
|
TEST_CASE("LatticeSearch - orthorhombic F - permutation 1") {
|
|
const double a = 41.0, b = 57.0, c = 35.0;
|
|
CrystalLattice conv(a,b,c, 90.0,90.0,90.0);
|
|
|
|
CrystalLattice L = conv.ToPrimitive('F');
|
|
|
|
auto res = LatticeSearch(L, 1e-6);
|
|
|
|
CHECK(res.system == gemmi::CrystalSystem::Orthorhombic);
|
|
CHECK(res.centering == 'F');
|
|
|
|
auto uc = res.conventional.GetUnitCell();
|
|
check_uc(uc, c, a, b, 90.0, 90.0, 90.0, 1e-1, 1e-2);
|
|
}
|
|
|
|
// Orthorhombic C: all angles 90, unequal edges, C-centering
|
|
TEST_CASE("LatticeSearch - orthorhombic C") {
|
|
const double a = 35.0, b = 41.0, c = 57.0;
|
|
CrystalLattice conv(a,b,c, 90.0,90.0,90.0);
|
|
CrystalLattice L = conv.ToPrimitive('C');
|
|
|
|
auto res = LatticeSearch(L, 1e-6);
|
|
|
|
CHECK(res.system == gemmi::CrystalSystem::Orthorhombic);
|
|
CHECK(res.centering == 'C');
|
|
|
|
auto uc = res.conventional.GetUnitCell();
|
|
check_uc(uc, a, b, c, 90.0, 90.0, 90.0, 1e-1, 1e-2);
|
|
}
|
|
|
|
TEST_CASE("LatticeSearch - orthorhombic I") {
|
|
const double a = 35.0, b = 41.0, c = 57.0;
|
|
CrystalLattice conv(a,b,c, 90.0,90.0,90.0);
|
|
CrystalLattice L = conv.ToPrimitive('I');
|
|
|
|
auto res = LatticeSearch(L, 1e-6);
|
|
|
|
CHECK(res.system == gemmi::CrystalSystem::Orthorhombic);
|
|
CHECK(res.centering == 'I');
|
|
|
|
auto uc = res.conventional.GetUnitCell();
|
|
check_uc(uc, a, b, c, 90.0, 90.0, 90.0, 1e-2, 1e-2);
|
|
}
|
|
|
|
|
|
TEST_CASE("LatticeSearch - orthorhombic I - permutation1") {
|
|
const double a = 57.0, b = 41.0, c = 35.0;
|
|
CrystalLattice conv(a,b,c, 90.0,90.0,90.0);
|
|
CrystalLattice L = conv.ToPrimitive('I');
|
|
|
|
auto res = LatticeSearch(L, 1e-6);
|
|
|
|
CHECK(res.system == gemmi::CrystalSystem::Orthorhombic);
|
|
CHECK(res.centering == 'I');
|
|
|
|
auto uc = res.conventional.GetUnitCell();
|
|
check_uc(uc, c, b, a, 90.0, 90.0, 90.0, 1e-2, 1e-2);
|
|
}
|
|
|
|
TEST_CASE("LatticeSearch - orthorhombic I - permutation2") {
|
|
const double a = 41.0, b = 57.0, c = 35.0;
|
|
CrystalLattice conv(a,b,c, 90.0,90.0,90.0);
|
|
CrystalLattice L = conv.ToPrimitive('I');
|
|
|
|
auto res = LatticeSearch(L, 1e-6);
|
|
|
|
CHECK(res.system == gemmi::CrystalSystem::Orthorhombic);
|
|
CHECK(res.centering == 'I');
|
|
|
|
auto uc = res.conventional.GetUnitCell();
|
|
check_uc(uc, c, a, b, 90.0, 90.0, 90.0, 1e-2, 1e-2);
|
|
}
|
|
|
|
// A character states its scalar products as fractions of A, B and C, and the three C-centred
|
|
// monoclinic ones (28, 29, 30) state one of them as 2*D or 2*E - twice a cosine. The cosine that
|
|
// implies leaves [-1,1] as soon as the cell's own angle is far enough from 90, and the character is
|
|
// then geometrically impossible for that metric. This cell is triclinic; character 28 asks it for a
|
|
// gamma whose cosine is 1.127.
|
|
TEST_CASE("LatticeSearch - an impossible character is not a match") {
|
|
CrystalLattice L(30.0, 35.0, 40.0, 65.0, 70.0, 70.0);
|
|
|
|
auto res = LatticeSearch(L, 1e-6);
|
|
|
|
CHECK(res.system == gemmi::CrystalSystem::Triclinic);
|
|
}
|
|
|
|
// An exact I-centred orthorhombic lattice whose reduced cell comes out all-acute with gamma at 90 -
|
|
// ON the boundary between the two Niggli types, where the reduction may present either. Character 42
|
|
// is stated for the obtuse setting, and only the flip that keeps gamma reaches it. Both defects have
|
|
// to be gone: without the impossible-character fix this metric matches character 28 and never gets
|
|
// as far as the retry, and without the gamma flip the retry does not have the setting it needs.
|
|
TEST_CASE("LatticeSearch - orthorhombic I on the type boundary in gamma") {
|
|
const double a = 45.0, b = 50.0, c = 80.0;
|
|
CrystalLattice conv(a, b, c, 90.0, 90.0, 90.0);
|
|
CrystalLattice L = conv.ToPrimitive('I');
|
|
|
|
auto res = LatticeSearch(L, 1e-6);
|
|
|
|
CHECK(res.system == gemmi::CrystalSystem::Orthorhombic);
|
|
CHECK(res.centering == 'I');
|
|
|
|
auto uc = res.conventional.GetUnitCell();
|
|
check_uc(uc, a, b, c, 90.0, 90.0, 90.0, 1e-2, 1e-2);
|
|
}
|
|
|
|
// Orthorhombic P: all angles 90, unequal edges, P-centering
|
|
TEST_CASE("LatticeSearch - orthorhombic P") {
|
|
const double a = 35.0, b = 41.0, c = 57.0;
|
|
CrystalLattice L(a,b,c, 90.0,90.0,90.0);
|
|
|
|
auto res = LatticeSearch(L, 1e-6);
|
|
|
|
CHECK(res.system == gemmi::CrystalSystem::Orthorhombic);
|
|
CHECK(res.centering == 'P');
|
|
|
|
auto uc = res.conventional.GetUnitCell();
|
|
check_uc(uc, a, b, c, 90.0, 90.0, 90.0, 1e-6, 1e-4);
|
|
}
|
|
|
|
// Hexagonal P: a=b!=c, alpha=beta=90, gamma=120, P-centering
|
|
TEST_CASE("LatticeSearch - hexagonal P") {
|
|
const double a = 30.0, c = 48.0;
|
|
CrystalLattice L(
|
|
Coord(a, 0, 0),
|
|
Coord(-a/2, a*std::sqrt(3)/2, 0),
|
|
Coord(0, 0, c)
|
|
);
|
|
|
|
auto res = LatticeSearch(L, 1e-6);
|
|
|
|
CHECK(res.system == gemmi::CrystalSystem::Hexagonal);
|
|
CHECK(res.centering == 'P');
|
|
|
|
auto uc = res.conventional.GetUnitCell();
|
|
CHECK(uc.a == Catch::Approx(a).margin(1e-2));
|
|
CHECK(uc.b == Catch::Approx(a).margin(1e-2));
|
|
CHECK(uc.c == Catch::Approx(c).margin(1e-2));
|
|
CHECK(uc.alpha == Catch::Approx(90.0).margin(1e-2));
|
|
CHECK(uc.beta == Catch::Approx(90.0).margin(1e-2));
|
|
CHECK(uc.gamma == Catch::Approx(120.0).margin(1e-2));
|
|
}
|
|
|
|
TEST_CASE("LatticeSearch - monoclinic C (unique b)") {
|
|
const double a = 50.0, b = 60.0, c = 70.0;
|
|
const double alpha = 90.0, beta = 96.0, gamma = 90.0;
|
|
CrystalLattice conv(a,b,c, alpha,beta,gamma);
|
|
auto L = conv.ToPrimitive('C');
|
|
|
|
auto res = LatticeSearch(L, 1e-6);
|
|
|
|
CHECK(res.system == gemmi::CrystalSystem::Monoclinic);
|
|
CHECK(res.centering == 'C');
|
|
|
|
auto uc = res.conventional.GetUnitCell();
|
|
// Check right angles at alpha,gamma and non-90 beta; lengths comparable
|
|
CHECK(std::fabs(uc.alpha - 90.0) < 1e-3);
|
|
CHECK(std::fabs(uc.gamma - 90.0) < 1e-3);
|
|
CHECK(std::fabs(uc.beta - beta) < 1e-2);
|
|
// Lengths should match within small tolerance
|
|
CHECK(uc.a == Catch::Approx(a).margin(1e-2));
|
|
CHECK(uc.b == Catch::Approx(b).margin(1e-2));
|
|
CHECK(uc.c == Catch::Approx(c).margin(1e-2));
|
|
}
|
|
|
|
TEST_CASE("LatticeSearch - monoclinic C (unique b) - v2") {
|
|
const double a = 71.0, b = 35.0, c = 90.0;
|
|
const double alpha = 90.0, beta = 96.0, gamma = 90.0;
|
|
CrystalLattice conv(a,b,c, alpha,beta,gamma);
|
|
auto L = conv.ToPrimitive('C');
|
|
|
|
auto res = LatticeSearch(L, 1e-6);
|
|
|
|
CHECK(res.system == gemmi::CrystalSystem::Monoclinic);
|
|
CHECK(res.centering == 'C');
|
|
|
|
auto uc = res.conventional.GetUnitCell();
|
|
// Check right angles at alpha,gamma and non-90 beta; lengths comparable
|
|
CHECK(std::fabs(uc.alpha - 90.0) < 1e-3);
|
|
CHECK(std::fabs(uc.gamma - 90.0) < 1e-3);
|
|
CHECK(std::fabs(uc.beta - beta) < 1e-2);
|
|
// Lengths should match within small tolerance
|
|
CHECK(uc.a == Catch::Approx(a).margin(1e-2));
|
|
CHECK(uc.b == Catch::Approx(b).margin(1e-2));
|
|
CHECK(uc.c == Catch::Approx(c).margin(1e-2));
|
|
}
|
|
|
|
TEST_CASE("LatticeSearch - monoclinic C (unique a)") {
|
|
const double a = 60.0, b = 50.0, c = 70.0;
|
|
const double alpha = 96.0, beta = 90.0, gamma = 90.0;
|
|
CrystalLattice conv(a,b,c, alpha,beta,gamma);
|
|
auto L = conv.ToPrimitive('C');
|
|
|
|
auto res = LatticeSearch(L, 1e-6);
|
|
|
|
CHECK(res.system == gemmi::CrystalSystem::Monoclinic);
|
|
CHECK(res.centering == 'C');
|
|
|
|
auto uc = res.conventional.GetUnitCell();
|
|
// Check right angles at alpha,gamma and non-90 beta; lengths comparable
|
|
CHECK(std::fabs(uc.alpha - 90.0) < 1e-3);
|
|
CHECK(std::fabs(uc.gamma - 90.0) < 1e-3);
|
|
CHECK(std::fabs(uc.beta - alpha) < 1e-2);
|
|
// Lengths should match within small tolerance
|
|
CHECK(uc.a == Catch::Approx(b).margin(1e-2));
|
|
CHECK(uc.b == Catch::Approx(a).margin(1e-2));
|
|
CHECK(uc.c == Catch::Approx(c).margin(1e-2));
|
|
}
|
|
|
|
TEST_CASE("LatticeSearch - monoclinic P (unique b)") {
|
|
const double a = 50.0, b = 60.0, c = 70.0;
|
|
const double alpha = 90.0, beta = 96.0, gamma = 90.0;
|
|
CrystalLattice conv(a,b,c, alpha,beta,gamma);
|
|
|
|
auto res = LatticeSearch(conv, 1e-6);
|
|
|
|
CHECK(res.system == gemmi::CrystalSystem::Monoclinic);
|
|
CHECK(res.centering == 'P');
|
|
|
|
auto uc = res.conventional.GetUnitCell();
|
|
// Check right angles at alpha,gamma and non-90 beta; lengths comparable
|
|
CHECK(std::fabs(uc.alpha - 90.0) < 1e-3);
|
|
CHECK(std::fabs(uc.gamma - 90.0) < 1e-3);
|
|
CHECK(std::fabs(uc.beta - beta) < 1e-2);
|
|
// Lengths should match within small tolerance
|
|
CHECK(uc.a == Catch::Approx(a).margin(1e-2));
|
|
CHECK(uc.b == Catch::Approx(b).margin(1e-2));
|
|
CHECK(uc.c == Catch::Approx(c).margin(1e-2));
|
|
}
|
|
|
|
TEST_CASE("LatticeSearch - monoclinic P (unique b) - v2") {
|
|
const double a = 90.0, b = 35.0, c = 71.0;
|
|
const double alpha = 90.0, beta = 96.0, gamma = 90.0;
|
|
CrystalLattice conv(a,b,c, alpha,beta,gamma);
|
|
|
|
auto res = LatticeSearch(conv, 1e-6);
|
|
|
|
CHECK(res.system == gemmi::CrystalSystem::Monoclinic);
|
|
CHECK(res.centering == 'P');
|
|
|
|
auto uc = res.conventional.GetUnitCell();
|
|
// Check right angles at alpha,gamma and non-90 beta; lengths comparable
|
|
CHECK(std::fabs(uc.alpha - 90.0) < 1e-3);
|
|
CHECK(std::fabs(uc.gamma - 90.0) < 1e-3);
|
|
CHECK(std::fabs(uc.beta - beta) < 1e-2);
|
|
// Lengths should match within small tolerance
|
|
CHECK(uc.a == Catch::Approx(c).margin(1e-2));
|
|
CHECK(uc.b == Catch::Approx(b).margin(1e-2));
|
|
CHECK(uc.c == Catch::Approx(a).margin(1e-2));
|
|
}
|
|
|
|
TEST_CASE("LatticeSearch - triclinic P") {
|
|
// General triclinic primitive cell
|
|
CrystalLattice L(33.1, 41.7, 52.3, 89.1, 85.0, 76.3);
|
|
|
|
auto res = LatticeSearch(L, 1e-6);
|
|
|
|
// System should be triclinic, centering P, and conventional equals some standardized primitive
|
|
CHECK(res.system == gemmi::CrystalSystem::Triclinic);
|
|
CHECK(res.centering == 'P');
|
|
|
|
// The conventional cell should be metric-equivalent to input. We verify only the system and centering here.
|
|
// Reduced primitive must be non-singular
|
|
auto uc_red = res.primitive_reduced.GetUnitCell();
|
|
CHECK(uc_red.a > 0);
|
|
CHECK(uc_red.b > 0);
|
|
CHECK(uc_red.c > 0);
|
|
}
|
|
|
|
TEST_CASE("LatticeSearch - triclinic P - v2") {
|
|
// General triclinic primitive cell
|
|
CrystalLattice L(33.1, 41.7, 52.3, 100, 92, 115);
|
|
|
|
auto res = LatticeSearch(L, 1e-6);
|
|
|
|
// System should be triclinic, centering P, and conventional equals some standardized primitive
|
|
CHECK(res.system == gemmi::CrystalSystem::Triclinic);
|
|
CHECK(res.centering == 'P');
|
|
|
|
// The conventional cell should be metric-equivalent to input. We verify only the system and centering here.
|
|
// Reduced primitive must be non-singular
|
|
auto uc_red = res.primitive_reduced.GetUnitCell();
|
|
CHECK(uc_red.a > 0);
|
|
CHECK(uc_red.b > 0);
|
|
CHECK(uc_red.c > 0);
|
|
}
|
|
|
|
TEST_CASE("LatticeSearch - trigonal R") {
|
|
const double a = 32.0;
|
|
const double alpha = 80.0;
|
|
|
|
// Build rhombohedral in rhombohedral setting (primitive axes a=b=c, alpha=beta=gamma)
|
|
CrystalLattice L(a, a, a, alpha, alpha, alpha);
|
|
|
|
auto res = LatticeSearch(L, 1e-6);
|
|
|
|
CHECK(res.system == gemmi::CrystalSystem::Trigonal);
|
|
CHECK(res.centering == 'R');
|
|
|
|
auto uc_red = res.conventional.GetUnitCell();
|
|
CHECK(uc_red.alpha == Catch::Approx(90).margin(1e-2));
|
|
CHECK(uc_red.beta == Catch::Approx(90).margin(1e-2));
|
|
CHECK(uc_red.gamma == Catch::Approx(120).margin(1e-2));
|
|
|
|
auto uc_prim = res.primitive_reduced.GetUnitCell();
|
|
CHECK(uc_prim.alpha == Catch::Approx(alpha).margin(1e-2));
|
|
CHECK(uc_prim.beta == Catch::Approx(alpha).margin(1e-2));
|
|
CHECK(uc_prim.gamma == Catch::Approx(alpha).margin(1e-2));
|
|
}
|
|
|
|
// The class-filtered walk: the same table, restricted to one Bravais class. A tetragonal-P lattice is
|
|
// also a C-centred orthorhombic one (a_C = a+b, b_C = -a+b, c_C = c), and asking for that class has to
|
|
// return that setting even though the plain search rightly prefers the tetragonal one.
|
|
TEST_CASE("LatticeSearchForClass - tetragonal P also has a C-centred orthorhombic setting") {
|
|
const double a = 50.0, c = 120.0;
|
|
const CrystalLattice L(a, a, c, 90, 90, 90);
|
|
|
|
const auto plain = LatticeSearch(L, 1e-6);
|
|
CHECK(plain.system == gemmi::CrystalSystem::Tetragonal);
|
|
CHECK(plain.centering == 'P');
|
|
|
|
const auto ortho = LatticeSearchForClass(L, gemmi::CrystalSystem::Orthorhombic, 'C', 1e-6);
|
|
REQUIRE(ortho.has_value());
|
|
CHECK(ortho->system == gemmi::CrystalSystem::Orthorhombic);
|
|
CHECK(ortho->centering == 'C');
|
|
const auto uc = ortho->conventional.GetUnitCell();
|
|
// The C cell is the face diagonal on a and b, so twice the volume and a = b = a_tet * sqrt(2).
|
|
CHECK(uc.a == Catch::Approx(a * std::sqrt(2.0)).margin(1e-4));
|
|
CHECK(uc.b == Catch::Approx(a * std::sqrt(2.0)).margin(1e-4));
|
|
CHECK(uc.c == Catch::Approx(c).margin(1e-4));
|
|
CHECK(uc.alpha == Catch::Approx(90).margin(1e-4));
|
|
CHECK(uc.beta == Catch::Approx(90).margin(1e-4));
|
|
CHECK(uc.gamma == Catch::Approx(90).margin(1e-4));
|
|
}
|
|
|
|
TEST_CASE("LatticeSearchForClass - a class the metric cannot carry is refused") {
|
|
// A general triclinic metric has no monoclinic-C setting, and an F-centred cubic lattice has no
|
|
// hexagonal-P one (its hexagonal description is R-centred).
|
|
const CrystalLattice tri(41.0, 47.0, 53.0, 71.0, 83.0, 97.0);
|
|
CHECK_FALSE(LatticeSearchForClass(tri, gemmi::CrystalSystem::Monoclinic, 'C').has_value());
|
|
|
|
const double a = 60.0;
|
|
const auto cubic_f = CrystalLattice(a, a, a, 90, 90, 90).ToPrimitive('F');
|
|
CHECK(LatticeSearch(cubic_f, 1e-6).centering == 'F');
|
|
CHECK_FALSE(LatticeSearchForClass(cubic_f, gemmi::CrystalSystem::Hexagonal, 'P').has_value());
|
|
// ... but its rhombohedral setting is there, which is what makes the refusal above a real answer
|
|
// rather than an artefact of the filter.
|
|
const auto rhomb = LatticeSearchForClass(cubic_f, gemmi::CrystalSystem::Trigonal, 'R');
|
|
REQUIRE(rhomb.has_value());
|
|
CHECK(rhomb->centering == 'R');
|
|
}
|
|
|
|
TEST_CASE("LatticeSearchForClass - asking for what the plain search found returns the same setting") {
|
|
const double a = 40.0;
|
|
const auto L = CrystalLattice(a, a, a, 90, 90, 90).ToPrimitive('I');
|
|
const auto plain = LatticeSearch(L, 1e-6);
|
|
const auto filtered = LatticeSearchForClass(L, plain.system, plain.centering, 1e-6);
|
|
REQUIRE(filtered.has_value());
|
|
CHECK(filtered->niggli_class == plain.niggli_class);
|
|
check_uc(filtered->conventional.GetUnitCell(), a, a, a, 90, 90, 90, 1e-4, 1e-4);
|
|
}
|
|
|
|
// The reduction epsilon. An exactly body-centred tetragonal lattice with c > a*sqrt(2) reduces to a
|
|
// character whose gamma is 90 EXACTLY, so the scalar product that decides the Niggli type is
|
|
// structurally zero and what a float lattice carries there is rounding. Axis-aligned that rounding
|
|
// happens to vanish - which is why the two tetragonal-I cases above pass - but every lattice the
|
|
// pipeline classifies is a refined, ROTATED one, and rotating this one about its own 4-fold is
|
|
// enough to lose the 4-fold on 38 of 60 rotations.
|
|
TEST_CASE("LatticeSearch - a body-centred tetragonal lattice keeps its 4-fold once it is rotated") {
|
|
CrystalLattice L(Coord(40, 0, 0), Coord(0, 40, 0), Coord(0, 0, 90));
|
|
L = L.ToPrimitive('I').Multiply(RotMatrix(0.3f, Coord(0, 0, 1)));
|
|
const auto res = LatticeSearch(L, 1e-6);
|
|
CHECK(res.system == gemmi::CrystalSystem::Tetragonal);
|
|
CHECK(res.centering == 'I');
|
|
}
|
|
|
|
// ITA character 43, the mI form. An ordinary centred-monoclinic crystal that happens to reduce into
|
|
// the form the table names mI - the same Bravais lattice in another setting, there is no fifteenth
|
|
// type. With that row absent the walk reaches character 44 and the centring is lost outright. The
|
|
// three monoclinic-C cases above reduce to characters 14, 39 and 14, so none of them samples it.
|
|
TEST_CASE("LatticeSearch - a centred monoclinic lattice that reduces to the mI form keeps its centring") {
|
|
const CrystalLattice L = CrystalLattice(35, 60, 30, 90, 120, 90).ToPrimitive('C');
|
|
const auto res = LatticeSearch(L, 1e-6);
|
|
CHECK(res.system == gemmi::CrystalSystem::Monoclinic);
|
|
CHECK(res.centering == 'I');
|
|
const auto uc = res.conventional.GetUnitCell();
|
|
CHECK(std::fabs(uc.alpha - 90.0) < 1e-3);
|
|
CHECK(std::fabs(uc.gamma - 90.0) < 1e-3);
|
|
CHECK(std::fabs(res.conventional.CalcVolume())
|
|
== Catch::Approx(2 * std::fabs(L.CalcVolume())).epsilon(1e-4));
|
|
}
|