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
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
|
// ************************************************************************************************
//
// BornAgain: simulate and fit reflection and scattering
//
//! @file Sample/Lattice/Lattice3D.cpp
//! @brief Implements class Lattice.
//!
//! @homepage http://www.bornagainproject.org
//! @license GNU General Public License v3 or higher (see COPYING)
//! @copyright Forschungszentrum Jülich GmbH 2018
//! @authors Scientific Computing Group at MLZ (see CITATION, AUTHORS)
//
// ************************************************************************************************
#include "Sample/Lattice/Lattice3D.h"
#include "Sample/Lattice/ISelectionRule.h"
#include <gsl/gsl_linalg.h>
#include <numbers>
using std::numbers::pi;
Lattice3D::Lattice3D(const R3& a, const R3& b, const R3& c)
: m_a(a)
, m_b(b)
, m_c(c)
{
computeReciprocalVectors();
}
Lattice3D::Lattice3D(const Lattice3D& lattice)
: Lattice3D(lattice.m_a, lattice.m_b, lattice.m_c)
{
if (lattice.m_selection_rule)
setSelectionRule(*lattice.m_selection_rule);
}
Lattice3D::~Lattice3D() = default;
Lattice3D Lattice3D::rotated(const RotMatrix& rotMatrix) const
{
R3 q1 = rotMatrix.transformed(m_a);
R3 q2 = rotMatrix.transformed(m_b);
R3 q3 = rotMatrix.transformed(m_c);
Lattice3D result = {q1, q2, q3};
if (m_selection_rule)
result.setSelectionRule(*m_selection_rule);
return result;
}
//! Currently unused but may be useful for checks
R3 Lattice3D::getMillerDirection(double h, double k, double l) const
{
R3 direction = h * m_ra + k * m_rb + l * m_rc;
return direction.unit_or_throw();
}
double Lattice3D::unitCellVolume() const
{
return std::abs(m_a.dot(m_b.cross(m_c)));
}
//! Currently only used in tests
void Lattice3D::reciprocalLatticeBasis(R3& ra, R3& rb, R3& rc) const
{
ra = m_ra;
rb = m_rb;
rc = m_rc;
}
I3 Lattice3D::nearestI3(const R3& q) const
{
return {(int)std::lround(q.dot(m_a) / (2 * pi)), (int)std::lround(q.dot(m_b) / (2 * pi)),
(int)std::lround(q.dot(m_c) / (2 * pi))};
}
std::vector<R3> Lattice3D::reciprocalLatticeVectorsWithinRadius(const R3& q, double dq) const
{
I3 nearest_coords = nearestI3(q);
int max_X = std::lround(m_a.mag() * dq / (2 * pi));
int max_Y = std::lround(m_b.mag() * dq / (2 * pi));
int max_Z = std::lround(m_c.mag() * dq / (2 * pi));
std::vector<R3> result;
for (int index_X = -max_X; index_X <= max_X; ++index_X) {
for (int index_Y = -max_Y; index_Y <= max_Y; ++index_Y) {
for (int index_Z = -max_Z; index_Z <= max_Z; ++index_Z) {
I3 coords = I3(index_X, index_Y, index_Z) + nearest_coords;
if (m_selection_rule && !m_selection_rule->coordinateSelected(coords))
continue;
R3 latticePoint = coords.x() * m_ra + coords.y() * m_rb + coords.z() * m_rc;
if ((latticePoint - q).mag() <= dq)
result.push_back(latticePoint);
}
}
}
return result;
}
void Lattice3D::computeReciprocalVectors() const
{
R3 q23 = m_b.cross(m_c);
R3 q31 = m_c.cross(m_a);
R3 q12 = m_a.cross(m_b);
m_ra = (2 * pi) / m_a.dot(q23) * q23;
m_rb = (2 * pi) / m_b.dot(q31) * q31;
m_rc = (2 * pi) / m_c.dot(q12) * q12;
}
void Lattice3D::setSelectionRule(const ISelectionRule& selection_rule)
{
m_selection_rule.reset(selection_rule.clone());
}
bool operator==(const Lattice3D& left, const Lattice3D& right)
{
if (left.basisVectorA() != right.basisVectorA())
return false;
if (left.basisVectorB() != right.basisVectorB())
return false;
if (left.basisVectorC() != right.basisVectorC())
return false;
bool same_selection = (!left.selectionRule() && !right.selectionRule())
|| (left.selectionRule() && right.selectionRule()
&& (left.selectionRule()->isEqualTo(*right.selectionRule())));
return same_selection;
}
|
