BornAgain  1.19.79
Open-source research software to simulate and fit neutron and x-ray reflectometry and grazing-incidence small-angle scattering
Lattice3D Class Reference

Description

A Bravais lattice, characterized by three basis vectors, and optionally an ISelectionRule.

Definition at line 30 of file Lattice3D.h.

Inheritance diagram for Lattice3D:
[legend]
Collaboration diagram for Lattice3D:
[legend]

Public Member Functions

 Lattice3D (const Lattice3D &lattice)
 
 Lattice3D (R3 a, R3 b, R3 c)
 
 ~Lattice3D () override
 
R3 basisVectorA () const
 Returns basis vector a. More...
 
R3 basisVectorB () const
 Returns basis vector b. More...
 
R3 basisVectorC () const
 Returns basis vector c. More...
 
void checkNodeArgs () const
 Raises exception if a parameter value is invalid. More...
 
std::string className () const final
 Returns the class name, to be hard-coded in each leaf class that inherits from INode. More...
 
R3 getMillerDirection (double h, double k, double l) const
 Returns normalized direction corresponding to the given Miller indices. More...
 
I3 nearestReciprocalLatticeVectorCoordinates (R3 q) const
 Returns the nearest reciprocal lattice point from a given vector. More...
 
virtual std::vector< const INode * > nodeChildren () const
 Returns all children. More...
 
std::vector< const INode * > nodeOffspring () const
 Returns all descendants. More...
 
Lattice3Doperator= (const Lattice3D &)=delete
 
virtual std::vector< ParaMetaparDefs () const
 Returns the parameter definitions, to be hard-coded in each leaf class. More...
 
void reciprocalLatticeBasis (R3 &ra, R3 &rb, R3 &rc) const
 Returns the reciprocal basis vectors. More...
 
std::vector< R3 > reciprocalLatticeVectorsWithinRadius (R3 q, double dq) const
 Returns a list of reciprocal lattice vectors within distance dq of a vector q. More...
 
Lattice3D rotated (const RotMatrix &rotMatrix) const
 Creates rotated lattice. More...
 
void setSelectionRule (const ISelectionRule &selection_rule)
 Sets a selection rule for the reciprocal vectors. More...
 
double unitCellVolume () const
 Returns the volume of the unit cell. More...
 

Protected Attributes

std::vector< double > m_P
 

Private Member Functions

void computeReciprocalVectors () const
 

Private Attributes

R3 m_a
 
R3 m_b
 
R3 m_c
 Basis vectors in real space. More...
 
R3 m_ra
 
R3 m_rb
 
R3 m_rc
 Cache of basis vectors in reciprocal space. More...
 
std::unique_ptr< ISelectionRulem_selection_rule
 

Constructor & Destructor Documentation

◆ Lattice3D() [1/2]

Lattice3D::Lattice3D ( R3  a,
R3  b,
R3  c 
)

Definition at line 21 of file Lattice3D.cpp.

22  : m_a(a)
23  , m_b(b)
24  , m_c(c)
25 {
27 }
void computeReciprocalVectors() const
Definition: Lattice3D.cpp:99
R3 m_c
Basis vectors in real space.
Definition: Lattice3D.h:72

References computeReciprocalVectors().

Here is the call graph for this function:

◆ Lattice3D() [2/2]

Lattice3D::Lattice3D ( const Lattice3D lattice)

Definition at line 29 of file Lattice3D.cpp.

30  : Lattice3D(lattice.m_a, lattice.m_b, lattice.m_c)
31 {
32  if (lattice.m_selection_rule)
34 }
void setSelectionRule(const ISelectionRule &selection_rule)
Sets a selection rule for the reciprocal vectors.
Definition: Lattice3D.cpp:109
std::unique_ptr< ISelectionRule > m_selection_rule
Definition: Lattice3D.h:73
Lattice3D(R3 a, R3 b, R3 c)
Definition: Lattice3D.cpp:21

References m_selection_rule, and setSelectionRule().

Here is the call graph for this function:

◆ ~Lattice3D()

Lattice3D::~Lattice3D ( )
overridedefault

Member Function Documentation

◆ basisVectorA()

R3 Lattice3D::basisVectorA ( ) const
inline

Returns basis vector a.

Definition at line 43 of file Lattice3D.h.

43 { return m_a; }

References m_a.

Referenced by Interference3DLattice::Interference3DLattice(), and ReMesocrystal::calculateLargestReciprocalDistance().

◆ basisVectorB()

R3 Lattice3D::basisVectorB ( ) const
inline

Returns basis vector b.

Definition at line 46 of file Lattice3D.h.

46 { return m_b; }

References m_b.

Referenced by Interference3DLattice::Interference3DLattice(), and ReMesocrystal::calculateLargestReciprocalDistance().

◆ basisVectorC()

R3 Lattice3D::basisVectorC ( ) const
inline

Returns basis vector c.

Definition at line 49 of file Lattice3D.h.

49 { return m_c; }

References m_c.

Referenced by Interference3DLattice::Interference3DLattice(), and ReMesocrystal::calculateLargestReciprocalDistance().

◆ checkNodeArgs()

void INode::checkNodeArgs ( ) const
inherited

Raises exception if a parameter value is invalid.

Definition at line 27 of file INode.cpp.

28 {
29  size_t nP = m_P.size();
30  if (parDefs().size() != nP) {
31  std::cerr << "BUG in class " << className() << std::endl;
32  std::cerr << "#m_P = " << nP << std::endl;
33  std::cerr << "#PDf = " << parDefs().size() << std::endl;
34  for (const ParaMeta& pm : parDefs())
35  std::cerr << " PDf: " << pm.name << std::endl;
36  ASSERT(0);
37  }
38  ASSERT(parDefs().size() == nP);
39  for (size_t i = 0; i < nP; ++i) {
40  const ParaMeta pm = parDefs()[i];
41 
43  if (pm.vMin == -INF) {
44  ASSERT(pm.vMax == +INF);
45  // nothing to do
46  } else if (pm.vMax == +INF) {
47  ASSERT(pm.vMin == 0);
48  limits = RealLimits::nonnegative();
49  } else {
50  limits = RealLimits::limited(pm.vMin, pm.vMax);
51  }
52  limits.check(pm.name, m_P[i]);
53  }
54 }
#define ASSERT(condition)
Definition: Assert.h:45
const double INF
Definition: INode.h:26
virtual std::vector< ParaMeta > parDefs() const
Returns the parameter definitions, to be hard-coded in each leaf class.
Definition: INode.h:51
std::vector< double > m_P
Definition: INode.h:63
virtual std::string className() const =0
Returns the class name, to be hard-coded in each leaf class that inherits from INode.
Limits for a real fit parameter.
Definition: RealLimits.h:24
static RealLimits limitless()
Creates an object without bounds (default)
Definition: RealLimits.cpp:139
void check(const std::string &name, double value) const
Throws if value is outside limits. Parameter 'name' is for exception message.
Definition: RealLimits.cpp:170
static RealLimits nonnegative()
Creates an object which can have only positive values with 0. included.
Definition: RealLimits.cpp:124
static RealLimits limited(double left_bound_value, double right_bound_value)
Creates an object bounded from the left and right.
Definition: RealLimits.cpp:134
Metadata of one model parameter.
Definition: INode.h:29
double vMin
Definition: INode.h:33
double vMax
Definition: INode.h:34
std::string name
Definition: INode.h:30

References ASSERT, RealLimits::check(), INode::className(), INF, RealLimits::limited(), RealLimits::limitless(), INode::m_P, ParaMeta::name, RealLimits::nonnegative(), INode::parDefs(), ParaMeta::vMax, and ParaMeta::vMin.

Referenced by BarGauss::BarGauss(), BarLorentz::BarLorentz(), Bipyramid4::Bipyramid4(), Box::Box(), CantellatedCube::CantellatedCube(), Cone::Cone(), ConstantBackground::ConstantBackground(), CosineRippleBox::CosineRippleBox(), CosineRippleGauss::CosineRippleGauss(), CosineRippleLorentz::CosineRippleLorentz(), Cylinder::Cylinder(), DistributionCosine::DistributionCosine(), DistributionGate::DistributionGate(), DistributionGaussian::DistributionGaussian(), DistributionLogNormal::DistributionLogNormal(), DistributionLorentz::DistributionLorentz(), DistributionTrapezoid::DistributionTrapezoid(), Dodecahedron::Dodecahedron(), EllipsoidalCylinder::EllipsoidalCylinder(), FootprintGauss::FootprintGauss(), FootprintSquare::FootprintSquare(), FuzzySphere::FuzzySphere(), GaussSphere::GaussSphere(), HemiEllipsoid::HemiEllipsoid(), HollowSphere::HollowSphere(), HorizontalCylinder::HorizontalCylinder(), Icosahedron::Icosahedron(), LongBoxGauss::LongBoxGauss(), LongBoxLorentz::LongBoxLorentz(), PlatonicOctahedron::PlatonicOctahedron(), PlatonicTetrahedron::PlatonicTetrahedron(), Prism3::Prism3(), Prism6::Prism6(), Profile1DCauchy::Profile1DCauchy(), Profile1DCosine::Profile1DCosine(), Profile1DGate::Profile1DGate(), Profile1DGauss::Profile1DGauss(), Profile1DTriangle::Profile1DTriangle(), Profile1DVoigt::Profile1DVoigt(), Profile2DCauchy::Profile2DCauchy(), Profile2DCone::Profile2DCone(), Profile2DGate::Profile2DGate(), Profile2DGauss::Profile2DGauss(), Profile2DVoigt::Profile2DVoigt(), Pyramid2::Pyramid2(), Pyramid3::Pyramid3(), Pyramid4::Pyramid4(), Pyramid6::Pyramid6(), RotationEuler::RotationEuler(), RotationX::RotationX(), RotationY::RotationY(), RotationZ::RotationZ(), SawtoothRippleBox::SawtoothRippleBox(), SawtoothRippleGauss::SawtoothRippleGauss(), SawtoothRippleLorentz::SawtoothRippleLorentz(), Sphere::Sphere(), Spheroid::Spheroid(), TruncatedCube::TruncatedCube(), TruncatedSphere::TruncatedSphere(), and TruncatedSpheroid::TruncatedSpheroid().

Here is the call graph for this function:

◆ className()

std::string Lattice3D::className ( ) const
inlinefinalvirtual

Returns the class name, to be hard-coded in each leaf class that inherits from INode.

Implements INode.

Definition at line 37 of file Lattice3D.h.

37 { return "Lattice3D"; }

◆ computeReciprocalVectors()

void Lattice3D::computeReciprocalVectors ( ) const
private

Definition at line 99 of file Lattice3D.cpp.

100 {
101  R3 q23 = m_b.cross(m_c);
102  R3 q31 = m_c.cross(m_a);
103  R3 q12 = m_a.cross(m_b);
104  m_ra = M_TWOPI / m_a.dot(q23) * q23;
105  m_rb = M_TWOPI / m_b.dot(q31) * q31;
106  m_rc = M_TWOPI / m_c.dot(q12) * q12;
107 }
#define M_TWOPI
Definition: Constants.h:54
R3 m_rc
Cache of basis vectors in reciprocal space.
Definition: Lattice3D.h:75

References m_a, m_b, m_c, m_ra, m_rb, m_rc, and M_TWOPI.

Referenced by Lattice3D().

◆ getMillerDirection()

R3 Lattice3D::getMillerDirection ( double  h,
double  k,
double  l 
) const

Returns normalized direction corresponding to the given Miller indices.

Currently unused but may be useful for checks.

Definition at line 50 of file Lattice3D.cpp.

51 {
52  R3 direction = h * m_ra + k * m_rb + l * m_rc;
53  return direction.unit();
54 }

References m_ra, m_rb, and m_rc.

◆ nearestReciprocalLatticeVectorCoordinates()

I3 Lattice3D::nearestReciprocalLatticeVectorCoordinates ( R3  q) const

Returns the nearest reciprocal lattice point from a given vector.

Definition at line 69 of file Lattice3D.cpp.

70 {
71  return {(int)std::lround(q.dot(m_a) / M_TWOPI), (int)std::lround(q.dot(m_b) / M_TWOPI),
72  (int)std::lround(q.dot(m_c) / M_TWOPI)};
73 }

References m_a, m_b, m_c, and M_TWOPI.

Referenced by reciprocalLatticeVectorsWithinRadius().

◆ nodeChildren()

◆ nodeOffspring()

std::vector< const INode * > INode::nodeOffspring ( ) const
inherited

Returns all descendants.

Definition at line 61 of file INode.cpp.

62 {
63  std::vector<const INode*> result;
64  result.push_back(this);
65  for (const auto* child : nodeChildren()) {
66  for (const auto* p : child->nodeOffspring())
67  result.push_back(p);
68  }
69  return result;
70 }
virtual std::vector< const INode * > nodeChildren() const
Returns all children.
Definition: INode.cpp:56

References INode::nodeChildren().

Here is the call graph for this function:

◆ operator=()

Lattice3D& Lattice3D::operator= ( const Lattice3D )
delete

◆ parDefs()

virtual std::vector<ParaMeta> INode::parDefs ( ) const
inlinevirtualinherited

Returns the parameter definitions, to be hard-coded in each leaf class.

Reimplemented in ConstantBackground, GaussSphere, FuzzySphere, RotationEuler, RotationZ, RotationY, RotationX, Crystal, Layer, HexagonalLattice2D, SquareLattice2D, BasicLattice2D, LayerRoughness, TruncatedSpheroid, TruncatedSphere, TruncatedCube, Spheroid, Sphere, SawtoothRippleLorentz, SawtoothRippleGauss, SawtoothRippleBox, Pyramid6, Pyramid4, Pyramid3, Pyramid2, Prism6, Prism3, PlatonicTetrahedron, PlatonicOctahedron, LongBoxLorentz, LongBoxGauss, Icosahedron, HorizontalCylinder, HollowSphere, HemiEllipsoid, EllipsoidalCylinder, Dodecahedron, Cylinder, CosineRippleLorentz, CosineRippleGauss, CosineRippleBox, Cone, CantellatedCube, Box, Bipyramid4, BarLorentz, BarGauss, Profile2DVoigt, Profile2DCone, Profile2DGate, Profile2DGauss, Profile2DCauchy, Profile1DVoigt, Profile1DCosine, Profile1DTriangle, Profile1DGate, Profile1DGauss, Profile1DCauchy, MisesGaussPeakShape, MisesFisherGaussPeakShape, LorentzFisherPeakShape, GaussFisherPeakShape, IsotropicLorentzPeakShape, IsotropicGaussPeakShape, ParticleLayout, InterferenceTwin, InterferenceRadialParaCrystal, InterferenceHardDisk, Interference2DSuperLattice, Interference2DParaCrystal, Interference1DLattice, DistributionTrapezoid, DistributionCosine, DistributionLogNormal, DistributionGaussian, DistributionLorentz, DistributionGate, ResolutionFunction2DGaussian, PolFilter, FootprintSquare, and FootprintGauss.

Definition at line 51 of file INode.h.

51 { return {}; }

Referenced by INode::checkNodeArgs(), and IFormFactor::pythonConstructor().

◆ reciprocalLatticeBasis()

void Lattice3D::reciprocalLatticeBasis ( R3 &  ra,
R3 &  rb,
R3 &  rc 
) const

Returns the reciprocal basis vectors.

Currently only used in tests.

Definition at line 62 of file Lattice3D.cpp.

63 {
64  ra = m_ra;
65  rb = m_rb;
66  rc = m_rc;
67 }

References m_ra, m_rb, and m_rc.

◆ reciprocalLatticeVectorsWithinRadius()

std::vector< R3 > Lattice3D::reciprocalLatticeVectorsWithinRadius ( R3  q,
double  dq 
) const

Returns a list of reciprocal lattice vectors within distance dq of a vector q.

Definition at line 75 of file Lattice3D.cpp.

76 {
77  I3 nearest_coords = nearestReciprocalLatticeVectorCoordinates(q);
78 
79  int max_X = std::lround(m_a.mag() * dq / M_TWOPI);
80  int max_Y = std::lround(m_b.mag() * dq / M_TWOPI);
81  int max_Z = std::lround(m_c.mag() * dq / M_TWOPI);
82 
83  std::vector<R3> result;
84  for (int index_X = -max_X; index_X <= max_X; ++index_X) {
85  for (int index_Y = -max_Y; index_Y <= max_Y; ++index_Y) {
86  for (int index_Z = -max_Z; index_Z <= max_Z; ++index_Z) {
87  I3 coords = I3(index_X, index_Y, index_Z) + nearest_coords;
88  if (m_selection_rule && !m_selection_rule->coordinateSelected(coords))
89  continue;
90  R3 latticePoint = coords.x() * m_ra + coords.y() * m_rb + coords.z() * m_rc;
91  if ((latticePoint - q).mag() <= dq)
92  result.push_back(latticePoint);
93  }
94  }
95  }
96  return result;
97 }
I3 nearestReciprocalLatticeVectorCoordinates(R3 q) const
Returns the nearest reciprocal lattice point from a given vector.
Definition: Lattice3D.cpp:69

References m_a, m_b, m_c, m_ra, m_rb, m_rc, m_selection_rule, M_TWOPI, and nearestReciprocalLatticeVectorCoordinates().

Referenced by Interference3DLattice::iff_without_dw(), ReMesocrystal::theFF(), and ReMesocrystal::thePolFF().

Here is the call graph for this function:

◆ rotated()

Lattice3D Lattice3D::rotated ( const RotMatrix rotMatrix) const

Creates rotated lattice.

Definition at line 38 of file Lattice3D.cpp.

39 {
40  R3 q1 = rotMatrix.transformed(m_a);
41  R3 q2 = rotMatrix.transformed(m_b);
42  R3 q3 = rotMatrix.transformed(m_c);
43  Lattice3D result = {q1, q2, q3};
44  if (m_selection_rule)
46  return result;
47 }
A Bravais lattice, characterized by three basis vectors, and optionally an ISelectionRule.
Definition: Lattice3D.h:30
T transformed(const T &v) const
Return transformed vector v.
Definition: RotMatrix.cpp:76

References m_a, m_b, m_c, m_selection_rule, setSelectionRule(), and RotMatrix::transformed().

Here is the call graph for this function:

◆ setSelectionRule()

void Lattice3D::setSelectionRule ( const ISelectionRule selection_rule)

Sets a selection rule for the reciprocal vectors.

Definition at line 109 of file Lattice3D.cpp.

110 {
111  m_selection_rule.reset(selection_rule.clone());
112 }
virtual ISelectionRule * clone() const =0

References ISelectionRule::clone(), and m_selection_rule.

Referenced by Lattice3D(), and rotated().

Here is the call graph for this function:

◆ unitCellVolume()

double Lattice3D::unitCellVolume ( ) const

Returns the volume of the unit cell.

Definition at line 56 of file Lattice3D.cpp.

57 {
58  return std::abs(m_a.dot(m_b.cross(m_c)));
59 }

References m_a, m_b, and m_c.

Referenced by Compute::Slicing::createParticleInSlice(), ReMesocrystal::theFF(), and ReMesocrystal::thePolFF().

Member Data Documentation

◆ m_a

◆ m_b

◆ m_c

R3 Lattice3D::m_c
private

◆ m_P

◆ m_ra

R3 Lattice3D::m_ra
mutableprivate

◆ m_rb

R3 Lattice3D::m_rb
private

◆ m_rc

R3 Lattice3D::m_rc
private

Cache of basis vectors in reciprocal space.

Definition at line 75 of file Lattice3D.h.

Referenced by computeReciprocalVectors(), getMillerDirection(), reciprocalLatticeBasis(), and reciprocalLatticeVectorsWithinRadius().

◆ m_selection_rule

std::unique_ptr<ISelectionRule> Lattice3D::m_selection_rule
private

The documentation for this class was generated from the following files: