Drake
Drake C++ Documentation
Loading...
Searching...
No Matches
MonomialBasisElement Class Reference

Detailed Description

MonomialBasisElement represents a monomial, a product of powers of variables with non-negative integer exponents.

Note that it doesn't not include the coefficient part of a monomial. So x, x³y, xy²z are all valid MonomialBasisElement instances, but 1+x or 2xy²z are not. TODO(hongkai.dai): deprecate Monomial class and replace Monomial class with MonomialBasisElement class. For more information regarding the motivation of this class, please see Drake github issue #13602 and #13803.

#include <drake/common/symbolic/monomial_basis_element.h>

Public Member Functions

 MonomialBasisElement ()
 Constructs a monomial equal to 1.
 MonomialBasisElement (std::nullptr_t)
 Constructs a default value.
 MonomialBasisElement (const std::map< Variable, int > &var_to_degree_map)
 Constructs a MonomialBasisElement from variable to degree map.
 MonomialBasisElement (const Expression &e)
 Converts an expression to a monomial if the expression is written as ∏ᵢpow(xᵢ, kᵢ), otherwise throws a runtime error.
 MonomialBasisElement (const Eigen::Ref< const VectorX< Variable > > &vars, const Eigen::Ref< const Eigen::VectorXi > &degrees)
 Constructs a Monomial from a vector of variables vars and their corresponding integer degrees degrees.
 MonomialBasisElement (const Variable &var)
 Constructs a monomial basis element with only one variable, and the degree is 1.
 MonomialBasisElement (const Variable &var, int degree)
 Constructs a monomial basis element with only one variable, and the degree of that variable is given by degree.
 ~MonomialBasisElement () override
std::pair< double, MonomialBasisElementEvaluatePartial (const Environment &env) const
 Partially evaluates using a given environment env.
MonomialBasisElementpow_in_place (int p)
 Returns this monomial raised to p.
bool operator< (const MonomialBasisElement &other) const
 Compares two MonomialBasisElement in lexicographic order.
std::map< MonomialBasisElement, doubleDifferentiate (const Variable &var) const
 Differentiates this MonomialBasisElement.
std::map< MonomialBasisElement, doubleIntegrate (const Variable &var) const
 Integrates this MonomialBasisElement on a variable.
void MergeBasisElementInPlace (const MonomialBasisElement &other)
 Merges this basis element with another basis element other by merging their var_to_degree_map.
std::map< ChebyshevBasisElement, doubleToChebyshevBasis () const
 Converts this monomial to Chebyshev polynomial basis.
template<typename BasisElement>
std::map< BasisElement, doubleToBasis () const
 Converts this monomial to a weighted sum of basis elements of type BasisElement.
Implements CopyConstructible, CopyAssignable, MoveConstructible, MoveAssignable
 MonomialBasisElement (const MonomialBasisElement &)=default
MonomialBasisElementoperator= (const MonomialBasisElement &)=default
 MonomialBasisElement (MonomialBasisElement &&)=default
MonomialBasisElementoperator= (MonomialBasisElement &&)=default
Public Member Functions inherited from PolynomialBasisElement
 PolynomialBasisElement ()=default
 Constructs a polynomial basis with empty var_to_degree map.
 PolynomialBasisElement (const std::map< Variable, int > &var_to_degree_map)
 Constructs a polynomial basis given the variable and the degree of that variable.
 PolynomialBasisElement (const Eigen::Ref< const VectorX< Variable > > &vars, const Eigen::Ref< const Eigen::VectorXi > &degrees)
 Constructs a polynomial basis, such that it contains the variable-to-degree map vars(i)→degrees(i).
virtual ~PolynomialBasisElement ()
const std::map< Variable, int > & var_to_degree_map () const
const std::map< Variable, int > & get_powers () const
 Returns variable to degree map.
int total_degree () const
 Returns the total degree of a polynomial basis.
int degree (const Variable &v) const
 Returns the degree of this PolynomialBasisElement in a variable v.
Variables GetVariables () const
double Evaluate (const Environment &env) const
 Evaluates under a given environment env.
bool operator== (const PolynomialBasisElement &other) const
bool operator!= (const PolynomialBasisElement &other) const
Expression ToExpression () const
 PolynomialBasisElement (const PolynomialBasisElement &)=default
PolynomialBasisElementoperator= (const PolynomialBasisElement &)=default
 PolynomialBasisElement (PolynomialBasisElement &&)=default
PolynomialBasisElementoperator= (PolynomialBasisElement &&)=default

Friends

template<class HashAlgorithm>
void hash_append (HashAlgorithm &hasher, const MonomialBasisElement &item) noexcept
 Implements the hash_append generic hashing concept.

Additional Inherited Members

Protected Member Functions inherited from PolynomialBasisElement
bool lexicographical_compare (const PolynomialBasisElement &other) const
 Compares two PolynomialBasisElement using lexicographical order.
virtual bool EqualTo (const PolynomialBasisElement &other) const
void DoEvaluatePartial (const Environment &e, double *coeff, std::map< Variable, int > *new_basis_element) const
int * get_mutable_total_degree ()
std::map< Variable, int > * get_mutable_var_to_degree_map ()
void DoMergeBasisElementInPlace (const PolynomialBasisElement &other)
 Merge this basis element with another basis element by merging their var_to_degree_map.

Constructor & Destructor Documentation

◆ MonomialBasisElement() [1/9]

◆ MonomialBasisElement() [2/9]

◆ MonomialBasisElement() [3/9]

Constructs a monomial equal to 1.

Namely the toal degree is zero.

◆ MonomialBasisElement() [4/9]

MonomialBasisElement ( std::nullptr_t )
explicit

Constructs a default value.

This overload is used by Eigen when EIGEN_INITIALIZE_MATRICES_BY_ZERO is enabled.

◆ MonomialBasisElement() [5/9]

MonomialBasisElement ( const std::map< Variable, int > & var_to_degree_map)
explicit

Constructs a MonomialBasisElement from variable to degree map.

◆ MonomialBasisElement() [6/9]

MonomialBasisElement ( const Expression & e)
explicit

Converts an expression to a monomial if the expression is written as ∏ᵢpow(xᵢ, kᵢ), otherwise throws a runtime error.

Precondition
is_polynomial(e) should be true.

◆ MonomialBasisElement() [7/9]

MonomialBasisElement ( const Eigen::Ref< const VectorX< Variable > > & vars,
const Eigen::Ref< const Eigen::VectorXi > & degrees )

Constructs a Monomial from a vector of variables vars and their corresponding integer degrees degrees.

For example, MonomialBasisElement([x, y, z], [2, 0, 1]) constructs a MonomialBasisElement x²z.

Precondition
The size of vars should be the same as the size of degrees.
Exceptions
std::exceptionif degrees includes a negative integer.

◆ MonomialBasisElement() [8/9]

MonomialBasisElement ( const Variable & var)
explicit

Constructs a monomial basis element with only one variable, and the degree is 1.

◆ MonomialBasisElement() [9/9]

MonomialBasisElement ( const Variable & var,
int degree )

Constructs a monomial basis element with only one variable, and the degree of that variable is given by degree.

◆ ~MonomialBasisElement()

~MonomialBasisElement ( )
override

Member Function Documentation

◆ Differentiate()

std::map< MonomialBasisElement, double > Differentiate ( const Variable & var) const
nodiscard

Differentiates this MonomialBasisElement.

Since dxⁿ/dx = nxⁿ⁻¹, we return the map from the MonomialBasisElement to its coefficient. So if this MonomialBasisElement is x³y², then differentiate with x will return (x²y² → 3) as dx³y²/dx = 3x²y² If var is not a variable in MonomialBasisElement, then returns an empty map.

◆ EvaluatePartial()

std::pair< double, MonomialBasisElement > EvaluatePartial ( const Environment & env) const
nodiscard

Partially evaluates using a given environment env.

The evaluation result is of type pair<double, MonomialBasisElement>. The first component (: double) represents the coefficient part while the second component represents the remaining parts of the MonomialBasisElement which was not evaluated.

Example 1. Evaluate with a fully-specified environment (x³*y²).EvaluatePartial({{x, 2}, {y, 3}}) = (2³ * 3² = 8 * 9 = 72, MonomialBasisElement{} = 1).

Example 2. Evaluate with a partial environment (x³*y²).EvaluatePartial({{x, 2}}) = (2³ = 8, y²).

◆ Integrate()

std::map< MonomialBasisElement, double > Integrate ( const Variable & var) const
nodiscard

Integrates this MonomialBasisElement on a variable.

Since ∫ xⁿ dx = 1 / (n+1) xⁿ⁺¹, we return the map from the MonomialBasisElement to its coefficient in the integration result. So if this MonomialBasisElement is x³y², then we return (x⁴y² → 1/4) as ∫ x³y²dx = 1/4 x⁴y². If var is not a variable in this MonomialBasisElement, for example ∫ x³y²dz = x³y²z, then we return (x³y²z → 1)

◆ MergeBasisElementInPlace()

void MergeBasisElementInPlace ( const MonomialBasisElement & other)

Merges this basis element with another basis element other by merging their var_to_degree_map.

This is equivalent to multiplying this monomial basis element in place with monomial basis element other.

◆ operator<()

bool operator< ( const MonomialBasisElement & other) const

Compares two MonomialBasisElement in lexicographic order.

◆ operator=() [1/2]

MonomialBasisElement & operator= ( const MonomialBasisElement & )
default

◆ operator=() [2/2]

MonomialBasisElement & operator= ( MonomialBasisElement && )
default

◆ pow_in_place()

MonomialBasisElement & pow_in_place ( int p)

Returns this monomial raised to p.

Exceptions
std::exceptionif p is negative.

◆ ToBasis()

template<typename BasisElement>
std::map< BasisElement, double > ToBasis ( ) const

Converts this monomial to a weighted sum of basis elements of type BasisElement.

We return the map from each BasisElement to its coefficient. For example, if BasisElement=ChebyshevBasisElement, then when this = x²y³, it returns {[T₂(x)T₃(y)⇒1/8], [T₂(x)T₁(y)⇒3/8], [T₀(x)T₃(y)⇒1/8], [T₀(x)T₁(y)⇒3/8]}.

Note
Currently we only support
Template Parameters
BasisElementbeing MonomialBasisElement and ChebyshevBasisElement.

◆ ToChebyshevBasis()

std::map< ChebyshevBasisElement, double > ToChebyshevBasis ( ) const
nodiscard

Converts this monomial to Chebyshev polynomial basis.

For example,

  • For x², it returns 0.5T₂(x) + 0.5T₀(x).
  • For x²y³, it returns 1/8T₂(x)T₃(y) + 3/8T₂(x)T₁(y) + 1/8T₀(x)T₃(y) + 3/8T₀(x)T₁(y).

We return the map from each ChebyshevBasisElement to its coefficient. For example, when this = x², it returns {[T₂(x)⇒0.5], [T₀(x)⇒0.5]}. When this = x²y³, it returns {[T₂(x)T₃(y)⇒1/8], [T₂(x)T₁(y)⇒3/8], [T₀(x)T₃(y)⇒1/8], [T₀(x)T₁(y)⇒3/8]}.

◆ hash_append

template<class HashAlgorithm>
void hash_append ( HashAlgorithm & hasher,
const MonomialBasisElement & item )
friend

Implements the hash_append generic hashing concept.


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