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

Description

A cube, with tetrahedral truncation of all corners.

Definition at line 23 of file TruncatedCube.h.

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

Public Member Functions

 TruncatedCube (double length, double removed_length)
 
 TruncatedCube (std::vector< double > P)
 
void assert_platonic () const
 Assertions for Platonic solid. More...
 
double bottomZ (const IRotation *rotation) const override
 
virtual bool canSliceAnalytically (const IRotation *rot) const
 Default implementation only allows rotations along z-axis. 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...
 
TruncatedCubeclone () const override
 Returns a clone of this ISampleNode object. More...
 
std::vector< const Material * > containedMaterials () const
 Returns set of unique materials contained in this ISampleNode. More...
 
complex_t formfactor_at_bottom (C3 q) const override
 
virtual SpinMatrix formfactor_pol (C3 q) const
 Returns scattering amplitude for complex scattering wavevector q=k_i-k_f in case of matrix interactions. Default implementation calls formfactor_at_bottom(q) and multiplies with the unit matrix. More...
 
bool isMagnetic () const
 Returns true if there is any magnetic material in this ISampleNode. More...
 
double length () const
 
virtual const Materialmaterial () const
 Returns nullptr, unless overwritten to return a specific material. More...
 
virtual std::vector< const INode * > nodeChildren () const
 Returns all children. More...
 
std::vector< const INode * > nodeOffspring () const
 Returns all descendants. More...
 
std::vector< ParaMetaparDefs () const final
 Returns the parameter definitions, to be hard-coded in each leaf class. More...
 
virtual std::string pythonConstructor () const
 Creates the Python constructor of this class (or derived classes) More...
 
double radialExtension () const override
 Returns the (approximate in some cases) radial size of the particle of this form factor's shape. This is used for SSCA calculations. More...
 
double removedLength () const
 
std::string shapeName () const
 
virtual complex_t theFF (const WavevectorInfo &wavevectors) const
 
virtual SpinMatrix thePolFF (const WavevectorInfo &wavevectors) const
 
double topZ (const IRotation *rotation) const override
 
virtual void transferToCPP ()
 Used for Python overriding of clone (see swig/tweaks.py) More...
 
double volume () const override
 

Protected Member Functions

void setPolyhedron (const ff::PolyhedralTopology &topology, double z_bottom, const std::vector< R3 > &vertices)
 Called by child classes to set faces and other internal variables. More...
 

Protected Attributes

std::vector< double > m_P
 
std::unique_ptr< IShape3Dm_shape3D
 IShape3D object, used to retrieve vertices (which may be approximate in the case of round shapes). For soft particles, this will be a hard mean shape. More...
 

Private Attributes

const double & m_length
 
const double & m_removed_length
 
std::vector< R3 > m_vertices
 
double m_z_bottom
 for topZ, bottomZ computation only More...
 
std::unique_ptr< ff::Polyhedron > pimpl
 

Static Private Attributes

static const ff::PolyhedralTopology topology
 

Constructor & Destructor Documentation

◆ TruncatedCube() [1/2]

TruncatedCube::TruncatedCube ( double  length,
double  removed_length 
)

Definition at line 58 of file TruncatedCube.cpp.

59  : TruncatedCube(std::vector<double>{length, removed_length})
60 {
61 }
double length() const
Definition: TruncatedCube.h:38
TruncatedCube(double length, double removed_length)

References length().

Referenced by clone().

Here is the call graph for this function:

◆ TruncatedCube() [2/2]

TruncatedCube::TruncatedCube ( std::vector< double >  P)

Definition at line 33 of file TruncatedCube.cpp.

35  , m_length(m_P[0])
36  , m_removed_length(m_P[1])
37 {
38  checkNodeArgs();
39  if (m_removed_length > 0.5 * m_length) {
40  std::ostringstream ostr;
41  ostr << "::TruncatedCube() -> Error in class initialization ";
42  ostr << "with parameters 'length':" << m_length;
43  ostr << " 'removed_length':" << m_removed_length << "\n\n";
44  ostr << "Check for removed_length <= 0.5*length failed.";
45  throw std::runtime_error(ostr.str());
46  }
47 
48  double a = m_length / 2;
49  double c = a - m_removed_length;
50 
52  {{-c, -a, -a}, {-a, -c, -a}, {-a, -a, -c}, {c, -a, -a}, {a, -c, -a}, {a, -a, -c},
53  {-c, a, -a}, {-a, c, -a}, {-a, a, -c}, {c, a, -a}, {a, c, -a}, {a, a, -c},
54  {-c, -a, a}, {-a, -c, a}, {-a, -a, c}, {c, -a, a}, {a, -c, a}, {a, -a, c},
55  {-c, a, a}, {-a, c, a}, {-a, a, c}, {c, a, a}, {a, c, a}, {a, a, c}});
56 }
void setPolyhedron(const ff::PolyhedralTopology &topology, double z_bottom, const std::vector< R3 > &vertices)
Called by child classes to set faces and other internal variables.
IFormFactorPolyhedron(const std::vector< double > &PValues)
The mathematics implemented here is described in full detail in a paper by Joachim Wuttke,...
void checkNodeArgs() const
Raises exception if a parameter value is invalid.
Definition: INode.cpp:27
std::vector< double > m_P
Definition: INode.h:63
const double & m_removed_length
Definition: TruncatedCube.h:44
static const ff::PolyhedralTopology topology
Definition: TruncatedCube.h:42
const double & m_length
Definition: TruncatedCube.h:43

References INode::checkNodeArgs(), m_length, m_removed_length, IFormFactorPolyhedron::setPolyhedron(), and topology.

Here is the call graph for this function:

Member Function Documentation

◆ assert_platonic()

void IFormFactorPolyhedron::assert_platonic ( ) const
inherited

Assertions for Platonic solid.

Definition at line 79 of file IFormFactorPolyhedron.cpp.

80 {
81  pimpl->assert_platonic();
82 }
std::unique_ptr< ff::Polyhedron > pimpl

References IFormFactorPolyhedron::pimpl.

Referenced by Dodecahedron::Dodecahedron(), Icosahedron::Icosahedron(), PlatonicOctahedron::PlatonicOctahedron(), and PlatonicTetrahedron::PlatonicTetrahedron().

◆ bottomZ()

double IFormFactorPolyhedron::bottomZ ( const IRotation rotation) const
overridevirtualinherited

Reimplemented from IFormFactor.

Definition at line 53 of file IFormFactorPolyhedron.cpp.

54 {
55  return PolyhedralUtil::BottomZ(m_vertices, rotation);
56 }
std::vector< R3 > m_vertices
double BottomZ(const std::vector< R3 > &vertices, const IRotation *rotation)
Calculates the z-coordinate of the lowest vertex after rotation.

References PolyhedralUtil::BottomZ(), and IFormFactorPolyhedron::m_vertices.

Here is the call graph for this function:

◆ canSliceAnalytically()

bool IFormFactor::canSliceAnalytically ( const IRotation rot) const
virtualinherited

Default implementation only allows rotations along z-axis.

Reimplemented in Sphere.

Definition at line 64 of file IFormFactor.cpp.

65 {
66  return !rotation || rotation->zInvariant();
67 }

References IRotation::zInvariant().

Here is the call graph for this function:

◆ 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
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(), TruncatedSphere::TruncatedSphere(), and TruncatedSpheroid::TruncatedSpheroid().

Here is the call graph for this function:

◆ className()

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

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

Implements INode.

Definition at line 30 of file TruncatedCube.h.

30 { return "TruncatedCube"; }

◆ clone()

TruncatedCube* TruncatedCube::clone ( ) const
inlineoverridevirtual

Returns a clone of this ISampleNode object.

Implements IFormFactor.

Definition at line 29 of file TruncatedCube.h.

29 { return new TruncatedCube(m_length, m_removed_length); }

References TruncatedCube(), m_length, and m_removed_length.

Here is the call graph for this function:

◆ containedMaterials()

std::vector< const Material * > ISampleNode::containedMaterials ( ) const
inherited

Returns set of unique materials contained in this ISampleNode.

Definition at line 25 of file ISampleNode.cpp.

26 {
27  std::vector<const Material*> result;
28  if (const Material* p_material = material())
29  result.push_back(p_material);
30  for (const auto* child : nodeChildren()) {
31  if (const auto* sample = dynamic_cast<const ISampleNode*>(child)) {
32  for (const Material* p_material : sample->containedMaterials())
33  result.push_back(p_material);
34  }
35  }
36  return result;
37 }
virtual std::vector< const INode * > nodeChildren() const
Returns all children.
Definition: INode.cpp:56
Abstract base class for sample components and properties related to scattering.
Definition: ISampleNode.h:27
virtual const Material * material() const
Returns nullptr, unless overwritten to return a specific material.
Definition: ISampleNode.h:36
A wrapper for underlying material implementation.
Definition: Material.h:35

References ISampleNode::material(), and INode::nodeChildren().

Referenced by SampleUtils::Multilayer::ContainsCompatibleMaterials(), SampleToPython::initLabels(), and ISampleNode::isMagnetic().

Here is the call graph for this function:

◆ formfactor_at_bottom()

complex_t IFormFactorPolyhedron::formfactor_at_bottom ( C3  q) const
overridevirtualinherited

Implements IFormFactor.

Definition at line 63 of file IFormFactorPolyhedron.cpp.

64 {
65  return exp_I(-m_z_bottom * q.z()) * pimpl->formfactor(q);
66 }
double m_z_bottom
for topZ, bottomZ computation only

References IFormFactorPolyhedron::m_z_bottom, and IFormFactorPolyhedron::pimpl.

◆ formfactor_pol()

SpinMatrix IFormFactor::formfactor_pol ( C3  q) const
virtualinherited

Returns scattering amplitude for complex scattering wavevector q=k_i-k_f in case of matrix interactions. Default implementation calls formfactor_at_bottom(q) and multiplies with the unit matrix.

Definition at line 78 of file IFormFactor.cpp.

79 {
81 }
virtual complex_t formfactor_at_bottom(C3 q) const =0
static SpinMatrix One()
Definition: SpinMatrix.cpp:36

References IFormFactor::formfactor_at_bottom(), and SpinMatrix::One().

Referenced by IFormFactor::thePolFF().

Here is the call graph for this function:

◆ isMagnetic()

bool ISampleNode::isMagnetic ( ) const
inherited

Returns true if there is any magnetic material in this ISampleNode.

Definition at line 39 of file ISampleNode.cpp.

40 {
41  const auto materials = containedMaterials();
42  return std::any_of(materials.cbegin(), materials.cend(),
43  [](const Material* mat) { return mat->isMagneticMaterial(); });
44 }
std::vector< const Material * > containedMaterials() const
Returns set of unique materials contained in this ISampleNode.
Definition: ISampleNode.cpp:25

References ISampleNode::containedMaterials().

Referenced by reSample::make().

Here is the call graph for this function:

◆ length()

double TruncatedCube::length ( ) const
inline

Definition at line 38 of file TruncatedCube.h.

38 { return m_length; }

References m_length.

Referenced by TruncatedCube().

◆ material()

virtual const Material* ISampleNode::material ( ) const
inlinevirtualinherited

Returns nullptr, unless overwritten to return a specific material.

Reimplemented in Particle, and Layer.

Definition at line 36 of file ISampleNode.h.

36 { return nullptr; }

Referenced by ISampleNode::containedMaterials().

◆ 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 }

References INode::nodeChildren().

Here is the call graph for this function:

◆ parDefs()

std::vector<ParaMeta> TruncatedCube::parDefs ( ) const
inlinefinalvirtual

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

Reimplemented from INode.

Definition at line 32 of file TruncatedCube.h.

33  {
34  return {{"Length", "nm", "untruncated edge length", 0, +INF, 0},
35  {"RemovedLength", "nm", "edge length removed from one corner", 0, +INF, 0}};
36  }

References INF.

◆ pythonConstructor()

std::string IFormFactor::pythonConstructor ( ) const
virtualinherited

Creates the Python constructor of this class (or derived classes)

Definition at line 69 of file IFormFactor.cpp.

70 {
71  std::vector<std::pair<double, std::string>> arguments;
72  for (size_t i = 0; i < parDefs().size(); i++)
73  arguments.emplace_back(m_P[i], parDefs()[i].unit);
74 
75  return Py::Fmt::printFunction(className(), arguments);
76 }
std::string printFunction(const std::string &name, const std::vector< std::pair< double, std::string >> &arguments)
Print a function in the form "<name>(<arguments>)". arguments will be processed by printArguments(),...
Definition: PyFmt.cpp:168

References INode::className(), INode::m_P, INode::parDefs(), and Py::Fmt::printFunction().

Here is the call graph for this function:

◆ radialExtension()

double IFormFactorPolyhedron::radialExtension ( ) const
overridevirtualinherited

Returns the (approximate in some cases) radial size of the particle of this form factor's shape. This is used for SSCA calculations.

Implements IFormFactor.

Definition at line 72 of file IFormFactorPolyhedron.cpp.

73 {
74  return pimpl->radius();
75 }

References IFormFactorPolyhedron::pimpl.

◆ removedLength()

double TruncatedCube::removedLength ( ) const
inline

Definition at line 39 of file TruncatedCube.h.

39 { return m_removed_length; }

References m_removed_length.

◆ setPolyhedron()

void IFormFactorPolyhedron::setPolyhedron ( const ff::PolyhedralTopology &  topology,
double  z_bottom,
const std::vector< R3 > &  vertices 
)
protectedinherited

Called by child classes to set faces and other internal variables.

Definition at line 41 of file IFormFactorPolyhedron.cpp.

43 {
44  m_z_bottom = z_bottom;
45 
46  m_vertices.clear();
47  for (const R3& vertex : vertices)
48  m_vertices.push_back(vertex - R3{0, 0, z_bottom});
49 
50  pimpl = std::make_unique<ff::Polyhedron>(topology, vertices);
51 }

References IFormFactorPolyhedron::m_vertices, IFormFactorPolyhedron::m_z_bottom, and IFormFactorPolyhedron::pimpl.

Referenced by Bipyramid4::Bipyramid4(), CantellatedCube::CantellatedCube(), Dodecahedron::Dodecahedron(), Icosahedron::Icosahedron(), PlatonicOctahedron::PlatonicOctahedron(), PlatonicTetrahedron::PlatonicTetrahedron(), Pyramid2::Pyramid2(), Pyramid3::Pyramid3(), Pyramid4::Pyramid4(), Pyramid6::Pyramid6(), and TruncatedCube().

◆ shapeName()

std::string IFormFactor::shapeName ( ) const
inherited

Definition at line 33 of file IFormFactor.cpp.

34 {
35  if (className().substr(0, 10) == "FormFactor")
36  return className().substr(10);
37  return className();
38 }

References INode::className().

Here is the call graph for this function:

◆ theFF()

complex_t IFormFactor::theFF ( const WavevectorInfo wavevectors) const
virtualinherited

Definition at line 40 of file IFormFactor.cpp.

41 {
42  return formfactor_at_bottom(wavevectors.getQ());
43 }
C3 getQ() const

References IFormFactor::formfactor_at_bottom(), and WavevectorInfo::getQ().

Here is the call graph for this function:

◆ thePolFF()

SpinMatrix IFormFactor::thePolFF ( const WavevectorInfo wavevectors) const
virtualinherited

Definition at line 45 of file IFormFactor.cpp.

46 {
47  return formfactor_pol(wavevectors.getQ());
48 }
virtual SpinMatrix formfactor_pol(C3 q) const
Returns scattering amplitude for complex scattering wavevector q=k_i-k_f in case of matrix interactio...
Definition: IFormFactor.cpp:78

References IFormFactor::formfactor_pol(), and WavevectorInfo::getQ().

Here is the call graph for this function:

◆ topZ()

double IFormFactorPolyhedron::topZ ( const IRotation rotation) const
overridevirtualinherited

Reimplemented from IFormFactor.

Definition at line 58 of file IFormFactorPolyhedron.cpp.

59 {
60  return PolyhedralUtil::TopZ(m_vertices, rotation);
61 }
double TopZ(const std::vector< R3 > &vertices, const IRotation *rotation)
Calculates the z-coordinate of the highest vertex after rotation.

References IFormFactorPolyhedron::m_vertices, and PolyhedralUtil::TopZ().

Here is the call graph for this function:

◆ transferToCPP()

virtual void ICloneable::transferToCPP ( )
inlinevirtualinherited

Used for Python overriding of clone (see swig/tweaks.py)

Definition at line 32 of file ICloneable.h.

◆ volume()

double IFormFactorPolyhedron::volume ( ) const
overridevirtualinherited

Reimplemented from IFormFactor.

Definition at line 68 of file IFormFactorPolyhedron.cpp.

69 {
70  return pimpl->volume();
71 }

References IFormFactorPolyhedron::pimpl.

Member Data Documentation

◆ m_length

const double& TruncatedCube::m_length
private

Definition at line 43 of file TruncatedCube.h.

Referenced by TruncatedCube(), clone(), and length().

◆ m_P

◆ m_removed_length

const double& TruncatedCube::m_removed_length
private

Definition at line 44 of file TruncatedCube.h.

Referenced by TruncatedCube(), clone(), and removedLength().

◆ m_shape3D

◆ m_vertices

std::vector<R3> IFormFactorPolyhedron::m_vertices
privateinherited

◆ m_z_bottom

double IFormFactorPolyhedron::m_z_bottom
privateinherited

for topZ, bottomZ computation only

Definition at line 55 of file IFormFactorPolyhedron.h.

Referenced by IFormFactorPolyhedron::formfactor_at_bottom(), and IFormFactorPolyhedron::setPolyhedron().

◆ pimpl

◆ topology

const ff::PolyhedralTopology TruncatedCube::topology
staticprivate
Initial value:
= {{{{0, 1, 7, 6, 9, 10, 4, 3}, true},
{{0, 2, 1}, false},
{{3, 4, 5}, false},
{{9, 11, 10}, false},
{{6, 7, 8}, false},
{{0, 3, 5, 17, 15, 12, 14, 2}, true},
{{4, 10, 11, 23, 22, 16, 17, 5}, true},
{{1, 2, 14, 13, 19, 20, 8, 7}, true},
{{6, 8, 20, 18, 21, 23, 11, 9}, true},
{{15, 17, 16}, false},
{{12, 13, 14}, false},
{{18, 20, 19}, false},
{{21, 22, 23}, false},
{{12, 15, 16, 22, 21, 18, 19, 13}, true}},
true}

Definition at line 42 of file TruncatedCube.h.

Referenced by TruncatedCube().


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