Molassembler  3.0.0
Molecule graph and conformer library
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
Scine::Molassembler::BondStereopermutator Class Reference

Handles specific relative arrangements of two atom stereopermutators joined by a bond. More...

#include <BondStereopermutator.h>

Data Structures

struct  FittingReferences
 
struct  Impl
 

Public Types

Public types
enum  Alignment { Alignment::Eclipsed, Alignment::Staggered, Alignment::EclipsedAndStaggered, Alignment::BetweenEclipsedAndStaggered }
 How dihedrals are aligned in the generation of stereopermutations. More...
 
enum  FittingMode { FittingMode::Thresholded, FittingMode::Nearest }
 Differentiates how viable assignments are chosen during fitting. More...
 
using SitePositions = Eigen::Matrix< double, 3, Eigen::Dynamic >
 Spatial centroids of atom stereopermutator site atoms.
 
using SitePositionsPair = std::pair< SitePositions, SitePositions >
 Spatial centroids of two atom stereopermutators' site atoms.
 

Public Member Functions

Special member functions
 BondStereopermutator (BondStereopermutator &&other) noexcept
 
BondStereopermutatoroperator= (BondStereopermutator &&other) noexcept
 
 BondStereopermutator (const BondStereopermutator &other)
 
BondStereopermutatoroperator= (const BondStereopermutator &other)
 
 ~BondStereopermutator ()
 
Constructors
 BondStereopermutator ()=delete
 Constructs a bond stereopermutator on two atom stereopermutators without checking whether stereopermutations are feasible or not. More...
 
 BondStereopermutator (const AtomStereopermutator &stereopermutatorA, const AtomStereopermutator &stereopermutatorB, const BondIndex &edge, Alignment alignment=Alignment::Eclipsed)
 Constructs a bond stereopermutator on two atom stereopermutators without checking whether stereopermutations are feasible or not. More...
 
 BondStereopermutator (const PrivateGraph &graph, const StereopermutatorList &stereopermutators, const BondIndex &edge, Alignment alignment=Alignment::Eclipsed)
 Constructs a bond stereopermutator on two atom stereopermutators, removing obviously infeasible stereopermutations. More...
 
Modification
void assign (boost::optional< unsigned > assignment)
 Changes the assignment of the stereopermutator. More...
 
void assignRandom (Random::Engine &engine)
 Assign the Stereopermutator at random. More...
 
void applyPermutation (const std::vector< AtomIndex > &permutation)
 Applies an atom index permutation. More...
 
void fit (const SitePositionsPair &sitePositions, std::pair< FittingReferences, FittingReferences > fittingReferences, FittingMode mode=FittingMode::Thresholded)
 Determines the assignment the permutator is in from positional information. More...
 
void fit (const AngstromPositions &angstromWrapper, std::pair< FittingReferences, FittingReferences > fittingReferences, FittingMode mode=FittingMode::Thresholded)
 
void propagateGraphChange (const AtomStereopermutator::PropagatedState &oldPermutator, const AtomStereopermutator &newPermutator, const PrivateGraph &inner, const StereopermutatorList &permutators)
 Propagates the bond stereocenter's state through a possible ranking change on one of its constituting atom stereopermutators. More...
 
void propagateVertexRemoval (AtomIndex removedIndex)
 Propagates invalidated atom indices through the removal of an index.
 
Information
Alignment alignment () const
 Returns alignment parameter this was constructed with. More...
 
boost::optional< unsigned > assigned () const
 Returns the permutation index within the set of possible permutations, if set. More...
 
const
Stereopermutations::Composite
composite () const
 Gives read-only access to the underlying Composite object. More...
 
std::pair< AtomIndex, AtomIndexcompositeAlignment () const
 The atom identifier alignment of the permutational composite.
 
double dihedral (const AtomStereopermutator &stereopermutatorA, SiteIndex siteIndexA, const AtomStereopermutator &stereopermutatorB, SiteIndex siteIndexB) const
 Angle between sites at stereopermutators in the current assignment. More...
 
bool hasSameCompositeOrientation (const BondStereopermutator &other) const
 Returns whether this stereopermutator has the same relative orientation as another stereopermutator. More...
 
boost::optional< unsigned > indexOfPermutation () const
 Returns the index of permutation. More...
 
unsigned numAssignments () const
 Returns the number of possible assignments. More...
 
unsigned numStereopermutations () const
 Returns the number of possible stereopermutations. More...
 
std::string info () const
 Returns an information string for diagnostic purposes. More...
 
std::string rankInfo () const
 Returns an information for ranking equality checking purposes. More...
 
BondIndex placement () const
 Returns which bond this stereopermutator is placed on in the molecule. More...
 
Operators
bool operator< (const BondStereopermutator &other) const
 Compares whether the underlying composite and the assignment are equivalent to those of another bond stereopermutator.
 
bool operator== (const BondStereopermutator &other) const
 Compares whether the underlying composite and the assignment are equivalent to those of another bond stereopermutator.
 
bool operator!= (const BondStereopermutator &other) const
 Inverts operator==.
 

Static Public Member Functions

Static functions
static std::vector< unsigned > notObviouslyInfeasibleStereopermutations (const PrivateGraph &graph, const AtomStereopermutator &stereopermutatorA, const AtomStereopermutator &stereopermutatorB, const Stereopermutations::Composite &composite)
 Determine which stereopermutations aren't obviously infeasible. More...
 

Static Public Attributes

Static members
static constexpr double assignmentAcceptanceParameter = 0.1
 Assignment acceptance is dependent on this parameter. More...
 

Private Attributes

std::unique_ptr< ImplpImpl_
 

Detailed Description

Handles specific relative arrangements of two atom stereopermutators joined by a bond.

This class exists to model rotational barriers along bonds joining an arbitrary pair of idealized shapes.

Member Enumeration Documentation

How dihedrals are aligned in the generation of stereopermutations.

Enumerator
Eclipsed 

At least two shape vertices eclipse one another along the axis.

Staggered 

At least one pair of substituents are staggered along the axis.

EclipsedAndStaggered 

Both eclipsed and staggered alignments.

BetweenEclipsedAndStaggered 

Offset exactly halfway between eclipsed and staggered alignments.

Differentiates how viable assignments are chosen during fitting.

Enumerator
Thresholded 

Positions must be close to the idealized assignment geometry.

Nearest 

Choose whichever assignment best represents the geometry directly.

Constructor & Destructor Documentation

Scine::Molassembler::BondStereopermutator::BondStereopermutator ( )
delete

Constructs a bond stereopermutator on two atom stereopermutators without checking whether stereopermutations are feasible or not.

Complexity \(O(S!)\) where \(S\) is the size of the larger involved shape

Scine::Molassembler::BondStereopermutator::BondStereopermutator ( const AtomStereopermutator stereopermutatorA,
const AtomStereopermutator stereopermutatorB,
const BondIndex edge,
Alignment  alignment = Alignment::Eclipsed 
)

Constructs a bond stereopermutator on two atom stereopermutators without checking whether stereopermutations are feasible or not.

Complexity \(O(S!)\) where \(S\) is the size of the larger involved shape

Scine::Molassembler::BondStereopermutator::BondStereopermutator ( const PrivateGraph graph,
const StereopermutatorList stereopermutators,
const BondIndex edge,
Alignment  alignment = Alignment::Eclipsed 
)

Constructs a bond stereopermutator on two atom stereopermutators, removing obviously infeasible stereopermutations.

Complexity \(O(S!)\) where \(S\) is the size of the larger involved shape

Member Function Documentation

Alignment Scine::Molassembler::BondStereopermutator::alignment ( ) const

Returns alignment parameter this was constructed with.

Complexity \(\Theta(1)\)

void Scine::Molassembler::BondStereopermutator::applyPermutation ( const std::vector< AtomIndex > &  permutation)

Applies an atom index permutation.

Complexity \(\Theta(1)\)

Parameters
permutationThe permutation to apply
void Scine::Molassembler::BondStereopermutator::assign ( boost::optional< unsigned >  assignment)

Changes the assignment of the stereopermutator.

Parameters
assignmentThe new assignment of the stereopermutator. May be boost::none, which sets the chiral state as indeterminate. Must be less than the number of assignments if not None.

Complexity \(\Theta(1)\)

boost::optional<unsigned> Scine::Molassembler::BondStereopermutator::assigned ( ) const

Returns the permutation index within the set of possible permutations, if set.

Complexity \(\Theta(1)\)

Returns
whether the stereopermutator is assigned or not, and if so, which assignment it is.
void Scine::Molassembler::BondStereopermutator::assignRandom ( Random::Engine engine)

Assign the Stereopermutator at random.

Complexity \(\Theta(1)\)

Note
If the stereocenter is already assigned, it is reassigned.
The state of the passed PRNG is advanced.
const Stereopermutations::Composite& Scine::Molassembler::BondStereopermutator::composite ( ) const

Gives read-only access to the underlying Composite object.

Complexity \(\Theta(1)\)

double Scine::Molassembler::BondStereopermutator::dihedral ( const AtomStereopermutator stereopermutatorA,
SiteIndex  siteIndexA,
const AtomStereopermutator stereopermutatorB,
SiteIndex  siteIndexB 
) const

Angle between sites at stereopermutators in the current assignment.

Complexity \(\Theta(1)\)

Parameters
stereopermutatorAOne constituting stereopermutator
siteIndexAA site index of stereopermutatorA
stereopermutatorBThe other constituting stereopermutator
siteIndexBA site index of stereopermutatorB
Exceptions
std::logic_errorIf the stereopermutator is unassigned or if no dihedral is found for the passed sites
void Scine::Molassembler::BondStereopermutator::fit ( const SitePositionsPair sitePositions,
std::pair< FittingReferences, FittingReferences fittingReferences,
FittingMode  mode = FittingMode::Thresholded 
)

Determines the assignment the permutator is in from positional information.

The assignment of this permutator is determined from three-dimensional positions using penalties to modeled dihedral angles. An assignment is considered matched if each dihedral that is part of the permutation is within assignmentAcceptanceDihedralThreshold.

Complexity \(\Theta(P)\) where \(P\) is the number of permutations

Parameters
sitePositionsThe site positions of each constuting atom stereopermutator, following the same ordering as compositeAlignment()
fittingReferencesPair of FittingReferences to constituting atom stereopermutators (order is irrelevant)
modeMode altering the assignment of stereopermutations depending on geometric closeness to the idealized minimum
void Scine::Molassembler::BondStereopermutator::fit ( const AngstromPositions angstromWrapper,
std::pair< FittingReferences, FittingReferences fittingReferences,
FittingMode  mode = FittingMode::Thresholded 
)

This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts.

bool Scine::Molassembler::BondStereopermutator::hasSameCompositeOrientation ( const BondStereopermutator other) const

Returns whether this stereopermutator has the same relative orientation as another stereopermutator.

Complexity \(\Theta(1)\)

boost::optional<unsigned> Scine::Molassembler::BondStereopermutator::indexOfPermutation ( ) const

Returns the index of permutation.

Complexity \(\Theta(1)\)

std::string Scine::Molassembler::BondStereopermutator::info ( ) const

Returns an information string for diagnostic purposes.

Complexity \(\Theta(1)\)

static std::vector<unsigned> Scine::Molassembler::BondStereopermutator::notObviouslyInfeasibleStereopermutations ( const PrivateGraph graph,
const AtomStereopermutator stereopermutatorA,
const AtomStereopermutator stereopermutatorB,
const Stereopermutations::Composite composite 
)
static

Determine which stereopermutations aren't obviously infeasible.

Complexity \(\Theta(P)\), where \(P\) is the number of permutations in the underlying composite

unsigned Scine::Molassembler::BondStereopermutator::numAssignments ( ) const

Returns the number of possible assignments.

Complexity \(\Theta(1)\)

unsigned Scine::Molassembler::BondStereopermutator::numStereopermutations ( ) const

Returns the number of possible stereopermutations.

Complexity \(\Theta(1)\)

BondIndex Scine::Molassembler::BondStereopermutator::placement ( ) const

Returns which bond this stereopermutator is placed on in the molecule.

Complexity \(\Theta(1)\)

void Scine::Molassembler::BondStereopermutator::propagateGraphChange ( const AtomStereopermutator::PropagatedState oldPermutator,
const AtomStereopermutator newPermutator,
const PrivateGraph inner,
const StereopermutatorList permutators 
)

Propagates the bond stereocenter's state through a possible ranking change on one of its constituting atom stereopermutators.

Parameters
oldPermutatorThe AtomStereopermutator that is changing at its state before a ranking change effects its current stereopermutation or substituent ranking
newPermutatorThe AtomStereopermutator that is changing after its internal state has been propagated through a ranking change
innerInner graph for context
permutatorsStereopermutator list for context

Complexity \(O(S!)\) where \(S\) is the size of the larger involved shape

std::string Scine::Molassembler::BondStereopermutator::rankInfo ( ) const

Returns an information for ranking equality checking purposes.

Complexity \(\Theta(1)\)

Field Documentation

constexpr double Scine::Molassembler::BondStereopermutator::assignmentAcceptanceParameter = 0.1
static

Assignment acceptance is dependent on this parameter.

An assignment is accepted if the sum of absolute angular deviations of all dihedrals involved in a stereopermutation divided by the number of dihedrals is less than this factor times two pi divided by the larger number of involved substituents at the ends of the bond.

For the classic E/Z double bond case, there are two involved substituents at both sides. The tolerance per dihedral is then 0.1 * 360° / 2 = 18°.


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