From 83cf43b4e949390156c3fbd4e1d123f87ee86cdb Mon Sep 17 00:00:00 2001 From: "Michael J. Sparapany" Date: Mon, 29 Apr 2024 12:52:01 -0600 Subject: [PATCH] v0.3.2 --- CMakeLists.txt | 2 +- Lielab/__release__.hpp | 2 +- Lielab/domain/CompositeAlgebra.hpp | 82 ++++- Lielab/domain/CompositeManifold.hpp | 90 ++++- Lielab/domain/liealgebras.hpp | 1 + Lielab/domain/liealgebras/glc.hpp | 411 +++++++++++++++++++++++ Lielab/domain/liegroups.hpp | 1 + Lielab/domain/liegroups/CN.hpp | 2 +- Lielab/domain/liegroups/GLC.hpp | 267 +++++++++++++++ Lielab/domain/lieiii.hpp | 32 +- Lielab/functions/cayley1.hpp | 15 + Lielab/topos/functions/Ad.hpp | 5 + Lielab/topos/functions/cayley1.hpp | 4 + Lielab/topos/functions/dexp.hpp | 5 + Lielab/topos/functions/exp.hpp | 4 + Lielab/topos/functions/log.hpp | 4 + python/cppLielab.cpp | 101 ++++++ python/tests/domain/tests/test_domain.py | 80 +++++ tests/test_domain.cpp | 59 ++++ 19 files changed, 1128 insertions(+), 39 deletions(-) create mode 100644 Lielab/domain/liealgebras/glc.hpp create mode 100644 Lielab/domain/liegroups/GLC.hpp diff --git a/CMakeLists.txt b/CMakeLists.txt index 9738f82..bed054a 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -1,7 +1,7 @@ cmake_minimum_required(VERSION 3.23) project("Lielab" - VERSION 0.3.1 + VERSION 0.3.2 LANGUAGES CXX) if (CMAKE_PROJECT_NAME STREQUAL PROJECT_NAME) diff --git a/Lielab/__release__.hpp b/Lielab/__release__.hpp index 4b2b847..e8a6d93 100644 --- a/Lielab/__release__.hpp +++ b/Lielab/__release__.hpp @@ -5,7 +5,7 @@ namespace Lielab { -const std::string VERSION = "0.3.1"; +const std::string VERSION = "0.3.2"; const std::string AUTHOR = "Michael J. Sparapany"; const std::string CONTACT = "mjspara@sandia.gov"; const std::string LOCATION = "https://github.com/sandialabs/Lielab"; diff --git a/Lielab/domain/CompositeAlgebra.hpp b/Lielab/domain/CompositeAlgebra.hpp index 2aa1bcb..8cabd70 100644 --- a/Lielab/domain/CompositeAlgebra.hpp +++ b/Lielab/domain/CompositeAlgebra.hpp @@ -30,16 +30,18 @@ namespace domain class CompositeAlgebra { public: - static constexpr size_t INDEX_cn = 0; - static constexpr size_t INDEX_gl = 1; - static constexpr size_t INDEX_rn = 2; - static constexpr size_t INDEX_se = 3; - static constexpr size_t INDEX_so = 4; - static constexpr size_t INDEX_sp = 5; - static constexpr size_t INDEX_su = 6; + static constexpr size_t INDEX_cn = 0; + static constexpr size_t INDEX_gl = 1; + static constexpr size_t INDEX_glc = 2; + static constexpr size_t INDEX_rn = 3; + static constexpr size_t INDEX_se = 4; + static constexpr size_t INDEX_so = 5; + static constexpr size_t INDEX_sp = 6; + static constexpr size_t INDEX_su = 7; typedef std::variant(space[ii]).get_dimension(); } + else if (ind == INDEX_glc) + { + dim += std::get(space[ii]).get_dimension(); + } else if (ind == INDEX_rn) { dim += std::get(space[ii]).get_dimension(); @@ -136,6 +142,10 @@ class CompositeAlgebra { out(ii) = static_cast(std::get(space[ii]).shape); } + else if (ind == INDEX_glc) + { + out(ii) = static_cast(std::get(space[ii]).shape); + } else if (ind == INDEX_rn) { out(ii) = static_cast(std::get(space[ii]).shape); @@ -176,6 +186,10 @@ class CompositeAlgebra { vectors[ii] = std::get(this->space[ii]).get_vector(); } + else if (ind == INDEX_glc) + { + vectors[ii] = std::get(this->space[ii]).get_vector(); + } else if (ind == INDEX_rn) { vectors[ii] = std::get(this->space[ii]).get_vector(); @@ -219,6 +233,12 @@ class CompositeAlgebra std::get(this->space[ii]).set_vector(vec(Eigen::seqN(jj, dim))); jj += dim; } + else if (ind == INDEX_glc) + { + const size_t dim = std::get(this->space[ii]).get_dimension(); + std::get(this->space[ii]).set_vector(vec(Eigen::seqN(jj, dim))); + jj += dim; + } else if (ind == INDEX_rn) { const size_t dim = std::get(this->space[ii]).get_dimension(); @@ -269,6 +289,10 @@ class CompositeAlgebra { out.space.push_back(std::get(this->space[ii]) + std::get(other.space[ii])); } + else if (ind == INDEX_glc) + { + out.space.push_back(std::get(this->space[ii]) + std::get(other.space[ii])); + } else if (ind == INDEX_rn) { out.space.push_back(std::get(this->space[ii]) + std::get(other.space[ii])); @@ -309,6 +333,10 @@ class CompositeAlgebra { std::get(this->space[ii]) += std::get(other.space[ii]); } + else if (ind == INDEX_glc) + { + std::get(this->space[ii]) += std::get(other.space[ii]); + } else if (ind == INDEX_rn) { std::get(this->space[ii]) += std::get(other.space[ii]); @@ -349,6 +377,10 @@ class CompositeAlgebra { out.space.push_back(std::get(this->space[ii]) - std::get(other.space[ii])); } + else if (ind == INDEX_glc) + { + out.space.push_back(std::get(this->space[ii]) - std::get(other.space[ii])); + } else if (ind == INDEX_rn) { out.space.push_back(std::get(this->space[ii]) - std::get(other.space[ii])); @@ -389,6 +421,10 @@ class CompositeAlgebra { std::get(this->space[ii]) -= std::get(other.space[ii]); } + else if (ind == INDEX_glc) + { + std::get(this->space[ii]) -= std::get(other.space[ii]); + } else if (ind == INDEX_rn) { std::get(this->space[ii]) -= std::get(other.space[ii]); @@ -429,6 +465,10 @@ class CompositeAlgebra { out.space.push_back(-std::get(this->space[ii])); } + else if (ind == INDEX_glc) + { + out.space.push_back(-std::get(this->space[ii])); + } else if (ind == INDEX_rn) { out.space.push_back(-std::get(this->space[ii])); @@ -469,6 +509,10 @@ class CompositeAlgebra { out.space.push_back(std::get(this->space[ii]) * other); } + else if (ind == INDEX_glc) + { + out.space.push_back(std::get(this->space[ii]) * other); + } else if (ind == INDEX_rn) { out.space.push_back(std::get(this->space[ii]) * other); @@ -507,6 +551,10 @@ class CompositeAlgebra { std::get(this->space[ii]) *= other; } + else if (ind == INDEX_glc) + { + std::get(this->space[ii]) *= other; + } else if (ind == INDEX_rn) { std::get(this->space[ii]) *= other; @@ -547,6 +595,10 @@ class CompositeAlgebra { out.space.push_back(std::get(this->space[ii]) * std::get(other.space[ii])); } + else if (ind == INDEX_glc) + { + out.space.push_back(std::get(this->space[ii]) * std::get(other.space[ii])); + } else if (ind == INDEX_rn) { out.space.push_back(std::get(this->space[ii]) * std::get(other.space[ii])); @@ -587,6 +639,10 @@ class CompositeAlgebra { std::get(this->space[ii]) *= std::get(other.space[ii]); } + if (ind == INDEX_glc) + { + std::get(this->space[ii]) *= std::get(other.space[ii]); + } else if (ind == INDEX_rn) { std::get(this->space[ii]) *= std::get(other.space[ii]); @@ -627,6 +683,10 @@ class CompositeAlgebra { out.space.push_back(std::get(this->space[ii]) / other); } + else if (ind == INDEX_glc) + { + out.space.push_back(std::get(this->space[ii]) / other); + } else if (ind == INDEX_rn) { out.space.push_back(std::get(this->space[ii]) / other); @@ -665,6 +725,10 @@ class CompositeAlgebra { std::get(this->space[ii]) /= other; } + else if (ind == INDEX_glc) + { + std::get(this->space[ii]) /= other; + } else if (ind == INDEX_rn) { std::get(this->space[ii]) /= other; @@ -707,6 +771,10 @@ class CompositeAlgebra { out += "gl(" + std::to_string(shapes(ii)) + ")"; } + else if (ind == INDEX_glc) + { + out += "glc(" + std::to_string(shapes(ii)) + ")"; + } else if (ind == INDEX_rn) { out += "rn(" + std::to_string(shapes(ii)) + ")"; diff --git a/Lielab/domain/CompositeManifold.hpp b/Lielab/domain/CompositeManifold.hpp index b676f34..3c34e68 100644 --- a/Lielab/domain/CompositeManifold.hpp +++ b/Lielab/domain/CompositeManifold.hpp @@ -31,23 +31,26 @@ class CompositeManifold { public: - static constexpr size_t INDEX_CN = 0; - static constexpr size_t INDEX_GL = 1; - static constexpr size_t INDEX_RN = 2; - static constexpr size_t INDEX_SE = 3; - static constexpr size_t INDEX_SO = 4; - static constexpr size_t INDEX_SP = 5; - static constexpr size_t INDEX_SU = 6; - static constexpr size_t INDEX_cn = 7; - static constexpr size_t INDEX_gl = 8; - static constexpr size_t INDEX_rn = 9; - static constexpr size_t INDEX_se = 10; - static constexpr size_t INDEX_so = 11; - static constexpr size_t INDEX_sp = 12; - static constexpr size_t INDEX_su = 13; + static constexpr size_t INDEX_CN = 0; + static constexpr size_t INDEX_GL = 1; + static constexpr size_t INDEX_GLC = 2; + static constexpr size_t INDEX_RN = 3; + static constexpr size_t INDEX_SE = 4; + static constexpr size_t INDEX_SO = 5; + static constexpr size_t INDEX_SP = 6; + static constexpr size_t INDEX_SU = 7; + static constexpr size_t INDEX_cn = 8; + static constexpr size_t INDEX_gl = 9; + static constexpr size_t INDEX_glc = 10; + static constexpr size_t INDEX_rn = 11; + static constexpr size_t INDEX_se = 12; + static constexpr size_t INDEX_so = 13; + static constexpr size_t INDEX_sp = 14; + static constexpr size_t INDEX_su = 15; typedef std::variant(std::get(space[ii]).shape); } + else if (ind == INDEX_GLC) + { + out(ii) = static_cast(std::get(space[ii]).shape); + } else if (ind == INDEX_RN) { out(ii) = static_cast(std::get(space[ii]).shape); @@ -139,6 +147,10 @@ class CompositeManifold { out(ii) = static_cast(std::get(space[ii]).shape); } + else if (ind == INDEX_glc) + { + out(ii) = static_cast(std::get(space[ii]).shape); + } else if (ind == INDEX_rn) { out(ii) = static_cast(std::get(space[ii]).shape); @@ -187,6 +199,10 @@ class CompositeManifold { serials.push_back(std::get(space[ii]).serialize()); } + else if (ind == INDEX_GLC) + { + serials.push_back(std::get(space[ii]).serialize()); + } else if (ind == INDEX_RN) { serials.push_back(std::get(space[ii]).serialize()); @@ -215,6 +231,10 @@ class CompositeManifold { serials.push_back(std::get(space[ii]).get_vector()); } + else if (ind == INDEX_glc) + { + serials.push_back(std::get(space[ii]).get_vector()); + } else if (ind == INDEX_rn) { serials.push_back(std::get(space[ii]).get_vector()); @@ -261,6 +281,11 @@ class CompositeManifold sz = std::get(M).get_size(); std::get(M).unserialize(vec(Eigen::seqN(jj, sz))); } + else if (ind == INDEX_GLC) + { + sz = std::get(M).get_size(); + std::get(M).unserialize(vec(Eigen::seqN(jj, sz))); + } else if (ind == INDEX_RN) { sz = std::get(M).get_size(); @@ -296,6 +321,11 @@ class CompositeManifold sz = std::get(M).get_dimension(); std::get(M).set_vector(vec(Eigen::seqN(jj, sz))); } + else if (ind == INDEX_glc) + { + sz = std::get(M).get_dimension(); + std::get(M).set_vector(vec(Eigen::seqN(jj, sz))); + } else if (ind == INDEX_rn) { sz = std::get(M).get_dimension(); @@ -339,6 +369,10 @@ class CompositeManifold { out.space.push_back(std::get(this->space[ii]) * std::get(other.space[ii])); } + else if (ind == INDEX_GLC) + { + out.space.push_back(std::get(this->space[ii]) * std::get(other.space[ii])); + } else if (ind == INDEX_RN) { out.space.push_back(std::get(this->space[ii]) * std::get(other.space[ii])); @@ -367,6 +401,10 @@ class CompositeManifold { out.space.push_back(std::get(this->space[ii]) + std::get(other.space[ii])); } + else if (ind == INDEX_glc) + { + out.space.push_back(std::get(this->space[ii]) + std::get(other.space[ii])); + } else if (ind == INDEX_rn) { out.space.push_back(std::get(this->space[ii]) + std::get(other.space[ii])); @@ -407,6 +445,10 @@ class CompositeManifold { std::get(this->space[ii]) *= std::get(other.space[ii]); } + else if (ind == INDEX_GLC) + { + std::get(this->space[ii]) *= std::get(other.space[ii]); + } else if (ind == INDEX_RN) { std::get(this->space[ii]) *= std::get(other.space[ii]); @@ -435,6 +477,10 @@ class CompositeManifold { std::get(this->space[ii]) += std::get(other.space[ii]); } + else if (ind == INDEX_glc) + { + std::get(this->space[ii]) += std::get(other.space[ii]); + } else if (ind == INDEX_rn) { std::get(this->space[ii]) += std::get(other.space[ii]); @@ -475,6 +521,10 @@ class CompositeManifold { out.space.push_back(std::get(this->space[ii]).inverse()); } + if (ind == INDEX_GLC) + { + out.space.push_back(std::get(this->space[ii]).inverse()); + } else if (ind == INDEX_RN) { out.space.push_back(std::get(this->space[ii]).inverse()); @@ -503,6 +553,10 @@ class CompositeManifold { out.space.push_back(-std::get(this->space[ii])); } + else if (ind == INDEX_glc) + { + out.space.push_back(-std::get(this->space[ii])); + } else if (ind == INDEX_rn) { out.space.push_back(-std::get(this->space[ii])); @@ -545,6 +599,10 @@ class CompositeManifold { out += "GL(" + std::to_string(shapes(ii)) + ")"; } + else if (ind == INDEX_GLC) + { + out += "GLC(" + std::to_string(shapes(ii)) + ")"; + } else if (ind == INDEX_RN) { out += "RN(" + std::to_string(shapes(ii)) + ")"; @@ -573,6 +631,10 @@ class CompositeManifold { out += "gl(" + std::to_string(shapes(ii)) + ")"; } + else if (ind == INDEX_glc) + { + out += "glc(" + std::to_string(shapes(ii)) + ")"; + } else if (ind == INDEX_rn) { out += "rn(" + std::to_string(shapes(ii)) + ")"; diff --git a/Lielab/domain/liealgebras.hpp b/Lielab/domain/liealgebras.hpp index 1940c2b..9439418 100644 --- a/Lielab/domain/liealgebras.hpp +++ b/Lielab/domain/liealgebras.hpp @@ -2,6 +2,7 @@ #include "liealgebras/cn.hpp" #include "liealgebras/gl.hpp" +#include "liealgebras/glc.hpp" #include "liealgebras/rn.hpp" #include "liealgebras/se.hpp" #include "liealgebras/so.hpp" diff --git a/Lielab/domain/liealgebras/glc.hpp b/Lielab/domain/liealgebras/glc.hpp new file mode 100644 index 0000000..3253037 --- /dev/null +++ b/Lielab/domain/liealgebras/glc.hpp @@ -0,0 +1,411 @@ +#ifndef _LIELAB_DOMAIN_LIEALGEBRAS_glc_HPP +#define _LIELAB_DOMAIN_LIEALGEBRAS_glc_HPP + +#include "../../abstract/abstract_all.hpp" + +#include +#include + +namespace Lielab +{ +namespace domain +{ + +using Lielab::abstract::Imaginary; + +class glc +{ + public: + + /*! Flag specifying whether or not the algebra is abelian. */ + static constexpr bool abelian = false; + size_t shape = 0; + + + /*! Internal data for the algebra. */ + Eigen::MatrixXcd _data; + + glc() : glc(0) + { + /*! \f{equation*}{() \rightarrow \mathfrak{glc} \f} + * + * Empty initialization function. Enables instantiation like: + * + * Lielab::domain::glc x, y, z; + * + */ + + } + + glc(const size_t shape) + { + /*! \f{equation*}{(\mathbb{Z}) \rightarrow \mathfrak{glc} \f} + * + * Constructor instantiating an \f$\mathfrak{glc}\f$ object. + * + * Enables instantiation like: + * + * Lielab::domain::glc x(3), y(4), z(5); + * + * @param[in] shape The shape of the data matrix. + */ + + this->_data = Eigen::MatrixXcd::Zero(shape,shape); + this->shape = shape; + } + + glc(std::initializer_list> other) + { + /*! \f{equation}{(\mathbb{C}^{n \times 1}) \rightarrow \mathfrak{glc} \f} + * + * Constructor instantiating a \f$\mathfrak{glc}\f$ object from either an + * \f$n \times 1\f$ imaginary vector. + * + * @param[in] other The object to instantiate from as an imaginary vector. + */ + + const Eigen::VectorXcd vector_imag = Eigen::VectorXcd{std::move(other)}; + const ptrdiff_t sz = vector_imag.size(); + Eigen::VectorXd vector_real = Eigen::VectorXd::Zero(2*sz); + ptrdiff_t kk = 0; + for (ptrdiff_t ii = 0; ii < sz; ii++) + { + vector_real(kk) = std::real(vector_imag(ii)); + vector_real(kk+1) = std::imag(vector_imag(ii)); + kk += 2; + } + + this->set_vector(vector_real); + } + + template + glc(const Eigen::MatrixBase & other) + { + /*! \f{eqnarray*}{(\mathbb{C}^{n \times n}) &\rightarrow& \mathfrak{glc} \\ (\mathbb{C}^{n \times 1}) &\rightarrow& \mathfrak{glc} \f} + * + * Constructor instantiating an \f$\mathfrak{glc}\f$ object from an + * \f$n \times n\f$ imaginary matrix. + * + * @param[in] other The object to instantiate from as an imaginary matrix. + */ + + if (other.rows() == other.cols()) + { + shape = other.rows(); + _data = other; + } + else + { + throw Errorx("Input data to glc malformed."); + } + } + + template + glc & operator=(const Eigen::MatrixBase & other) + { + /*! \f{eqnarray*}{ \mathfrak{glc} &:= \mathbb{C}^{n \times n}\f} + * + * Overload of the assignment operator. Allows Eigen Matrix data to be directly assigned to `glc`. + */ + + if (other.rows() == other.cols()) + { + shape = other.rows(); + _data = other; + } + else + { + throw Errorx("Input data to u malformed."); + } + return *this; + } + + static glc basis(const size_t i, const size_t shape) + { + /*! \f{equation*}{ (\mathbb{Z}, \mathbb{Z}) \rightarrow \mathfrak{glc} \f} + * + * Returns the i'th basis element of the glc algebra. + * + * @param[in] i The basis vector. + * @param[in] shape The shape of the algebra. + * @param[out] out The glc element. + */ + + glc out(shape); + const size_t _dim = out.get_dimension(); + Eigen::VectorXd _u(_dim); + _u = Eigen::VectorXd::Zero(_dim); + _u(i) = 1.0; + out.set_vector(_u); + return out; + } + + static Eigen::MatrixXcd project(const Eigen::MatrixXcd & other) + { + /*! \f{equation*}{ (\mathbb{C}^{n \times n}) \rightarrow \mathbb{C}^{n \times n} \in \mathfrak{glc} \f} + * + * Projects a matrix suitable for data. + */ + + if (other.rows() != other.cols()) + { + throw Errorx("Size of the matrix must be square."); + } + return other; + } + + size_t get_dimension() const + { + /*! \f{equation*}{ () \rightarrow \mathbb{Z} \f} + * + * Gets the dimension of the algebra. + */ + + return static_cast(2*std::pow(this->shape, 2)); + } + + Eigen::VectorXd get_vector() const + { + /*! \f{equation*}{ () \rightarrow \mathbb{R}^{n \times 1} \f} + * + * Returns the vector representation. + */ + + const ptrdiff_t dim = this->get_dimension(); + const Eigen::MatrixXcd A = this->get_ados_representation(); + + Eigen::VectorXd out = Eigen::VectorXd::Zero(dim); + ptrdiff_t kk = 0; + for (ptrdiff_t ii = 0; ii < this->shape; ii++) + { + for (ptrdiff_t jj = 0; jj < this->shape; jj++) + { + out(kk) = std::real(A(ii,jj)); + out(kk+1) = std::imag(A(ii,jj)); + kk = kk + 2; + } + } + + return out; + } + + Eigen::MatrixXcd get_ados_representation() const + { + /*! \f{equation*}{ () \rightarrow \mathbb{C}^{n \times n} \f} + * + * Returns a matrix representation. + */ + + return _data; + } + + void set_vector(const Eigen::VectorXd & vector) + { + /*! \f{equation*}{ \mathfrak{glc} := \mathbb{R}^{n \times 1} \f} + * + * @param[in] vector An Eigen::VectorXd to assign. + */ + + const ptrdiff_t dim = vector.size(); + this->shape = static_cast(std::sqrt(dim/2)); + this->_data = Eigen::MatrixXcd::Zero(this->shape, this->shape); + + ptrdiff_t kk = 0; + for (ptrdiff_t ii = 0; ii < this->shape; ii++) + { + for (ptrdiff_t jj = 0; jj < this->shape; jj++) + { + this->_data(ii, jj) = std::complex(vector(kk), vector(kk+1)); + kk = kk + 2; + } + } + } + + double operator()(const size_t index) + { + /*! \f{equation*}{ (\mathbb{Z}) \rightarrow \mathbb{R} \f} + * + * Gets a value in the column vector representation. + */ + + const Eigen::VectorXd vector = get_vector(); + return vector(index); + } + + // double & operator()(const size_t index) // TODO [not important] + // { + // /*! \f{equation*}{ \mathfrak{glc}(\mathbb{Z}) := \mathbb{R} \f} + // * + // * Assignment of a value in the column vector representation. + // */ + + // return _data(index); + // } + + std::complex operator()(const size_t index1, const size_t index2) const + { + /*! \f{equation*}{ (\mathbb{Z}, \mathbb{Z}) \rightarrow \mathbb{C} \f} + * + * Gets a value in the square matrix representation. + */ + + return _data(index1, index2); + } + + std::complex & operator()(const size_t index1, const size_t index2) + { + /*! \f{equation*}{ \mathfrak{glc}(\mathbb{Z}, \mathbb{Z}) := \mathbb{C} \f} + * + * Assignment of a value in the square matrix representation. + */ + + return _data(index1, index2); + } + + glc operator+(const glc & other) const + { + /*! \f{equation*}{ (\mathfrak{glc}, \mathfrak{glc}) \rightarrow \mathfrak{glc} \f} + * + * Addition of two vectors in the algebra. + */ + + assert(this->shape == other.shape); + return this->_data + other._data; + } + + glc & operator+=(const glc & other) + { + /*! \f{equation*}{ (\mathfrak{glc}, \mathfrak{glc}) \rightarrow \mathfrak{glc} \f} + * + * In place addition of two vectors in the algebra. + */ + + assert(this->shape == other.shape); + this->_data += other._data; + return *this; + } + + glc operator-(const glc & other) const + { + /*! \f{equation*}{ (\mathfrak{glc}}, \mathfrak{glc}) \rightarrow \mathfrak{glc}} \f} + * + * Subtraction of two vectors in the algebra. + */ + + assert(this->shape == other.shape); + return this->_data - other._data; + } + + glc & operator-=(const glc & other) + { + /*! \f{equation*}{ (\mathfrak{glc}, \mathfrak{glc}) \rightarrow \mathfrak{glc} \f} + * + * In place subtraction of two vectors in the algebra. + */ + + assert(this->shape == other.shape); + this->_data -= other._data; + return *this; + } + + glc operator-() const + { + /*! \f{equation*}{ (\mathfrak{glc}) \rightarrow \mathfrak{glc} \f} + * + * Unary negative of the vector. + */ + + return -this->_data; + } + + glc operator*(const Imaginary auto other) const + { + /*! \f{equation*}{ (\mathfrak{glc}, \mathbb{C}) \rightarrow \mathfrak{glc} \f} + * + * Scalar product. + */ + + return this->_data * other; + } + + glc & operator*=(const Imaginary auto other) + { + /*! \f{equation*}{ (\mathfrak{glc}, \mathbb{C}) \rightarrow \mathfrak{glc} \f} + * + * In place scalar product. + */ + + this->_data *= other; + return *this; + } + + glc operator*(const glc & other) const + { + /*! \f{equation*}{ (\mathfrak{glc}, \mathfrak{glc}) \rightarrow \mathfrak{glc} \f} + * + * Vector product. + */ + + assert(this->shape == other.shape); + return this->_data * other._data; + } + + glc & operator*=(const glc & other) + { + /*! \f{equation*}{ (\mathfrak{glc}, \mathfrak{glc}) \rightarrow \mathfrak{glc} \f} + * + * In place product of two vectors in the algebra. + */ + + assert(this->shape == other.shape); + this->_data *= other._data; + return *this; + } + + glc operator/(const Imaginary auto other) const + { + /*! \f{equation*}{ (\mathfrak{glc}, \mathbb{C}) \rightarrow \mathfrak{glc} \f} + * + * Scalar division. + */ + + return this->_data / other; + } + + glc & operator/=(const Imaginary auto other) + { + /*! \f{equation*}{ (\mathfrak{glc}, \mathbb{C}) \rightarrow \mathfrak{glc} \f} + * + * In place scalar division. + */ + + this->_data /= other; + return *this; + } + + friend std::ostream& operator<<(std::ostream & os, const glc & other); +}; + +glc operator*(const Imaginary auto other, const glc & rhs) +{ + /*! \f{equation*}{ (\mathbb{C}, \mathfrak{glc}) \rightarrow \mathfrak{glc} \f} + * + * Scalar product. + */ + + return rhs*other; +} + +std::ostream& operator<<(std::ostream & os, const glc & other) +{ + /*! + * Overloads the "<<" stream insertion operator. + */ + + os << static_cast(other._data); + return os; +} +} +} + +#endif diff --git a/Lielab/domain/liegroups.hpp b/Lielab/domain/liegroups.hpp index 835be7e..1110deb 100644 --- a/Lielab/domain/liegroups.hpp +++ b/Lielab/domain/liegroups.hpp @@ -2,6 +2,7 @@ #include "liegroups/CN.hpp" #include "liegroups/GL.hpp" +#include "liegroups/GLC.hpp" #include "liegroups/RN.hpp" #include "liegroups/SE.hpp" #include "liegroups/SO.hpp" diff --git a/Lielab/domain/liegroups/CN.hpp b/Lielab/domain/liegroups/CN.hpp index e56312c..77481c4 100644 --- a/Lielab/domain/liegroups/CN.hpp +++ b/Lielab/domain/liegroups/CN.hpp @@ -170,7 +170,7 @@ namespace Lielab * Gets the size of the data representation. */ - return static_cast(this->shape - 1); + return static_cast(2*(this->shape - 1)); } Eigen::MatrixXcd get_ados_representation() const diff --git a/Lielab/domain/liegroups/GLC.hpp b/Lielab/domain/liegroups/GLC.hpp new file mode 100644 index 0000000..253150d --- /dev/null +++ b/Lielab/domain/liegroups/GLC.hpp @@ -0,0 +1,267 @@ +#ifndef _LIELAB_DOMAIN_GLC_HPP +#define _LIELAB_DOMAIN_GLC_HPP + +#include +#include + +namespace Lielab +{ + namespace domain + { + class GLC + { + /*! + * The GLC class. + */ + public: + static constexpr bool abelian = false; + size_t shape = 0; + + typedef Eigen::MatrixXcd Base; + Base _data; + + // Initialization methods + + GLC() : GLC(0) + { + /*! \f{equation*}{ () \rightarrow GLC \f} + * + * Empty initialization function. Enables instantiation like: + * + * Lielab::domain::GLC x, y, z; + * + */ + + } + + GLC(const size_t shape) + { + /*! \f{equation*}{ (\mathbb{Z}) \rightarrow GLC \f} + * + * Constructor instantiating an \f$GLC\f$ object. + * + * Enables instantiation like: + * + * Lielab::domain::GLC x(2), y(3), z(4); + * + * @param[in] shape The shape of the data matrix. + */ + + this->_data = Base::Identity(shape, shape); + this->shape = shape; + } + + template + GLC(const Eigen::MatrixBase & other) + { + /*! \f{eqnarray*}{(\mathbb{R}^{n \times n}) &\rightarrow& GLC \\ (\mathbb{R}^{n \times 1}) &\rightarrow& GLC \f} + * + * Constructor instantiating an \f$GLC\f$ object from either an + * \f$n \times n\f$ imaginary matrix or \f$n \times 1\f$ imaginary vector. + * + * @param[in] other The object to instantiate from as a real matrix. + */ + + if (other.rows() == other.cols()) + { + this->_data = other; + this->shape = other.rows(); + } + else + { + throw Errorx("Input data to GLC malformed."); + } + } + + template + GLC & operator=(const Eigen::MatrixBase & other) + { + /*! \f{eqnarray*}{ GLC &:= \mathbb{C}^{n \times n} \\ &:= \mathbb{C}^{n \times 1} \f} + * + * Overload of the assignment operator. Allows Eigen Matrix data to be directly assigned to `GLC`. + */ + + if (other.rows() == other.cols()) + { + shape = other.rows(); + _data = other; + } + return *this; + } + + static Eigen::MatrixXcd project(const Eigen::MatrixXcd & other) + { + /*! \f{equation*}{ (\mathbb{C}^{n \times n}) \rightarrow \mathbb{C}^{n \times n} \in GLC \f} + * + * Projects a matrix suitable for data. + */ + + if (other.rows() != other.cols()) + { + throw Errorx("Size of the matrix must be square."); + } + + return other; + } + + size_t get_dimension() const + { + /*! \f{equation*}{ () \rightarrow \mathbb{Z} \f} + * + * Gets the dimension of the group. + */ + + return static_cast(2*std::pow(this->shape, 2)); + } + + size_t get_size() const + { + /*! \f{quation*}{ () \rightarrow \mathbb{Z} \f} + * + * Gets the size of the data representation. + */ + + return static_cast(2*std::pow(this->shape, 2)); + } + + Eigen::MatrixXcd get_ados_representation() const + { + /*! \f{equation*}{ () \rightarrow \mathbb{C}^{n \times n} \f} + * + * Returns a matrix representation. + */ + + return _data; + } + + GLC inverse() const + { + /*! \f{equation*}{ (GLC) \rightarrow GLC \f} + * + * Returns the inverse. + */ + + return this->_data.inverse(); + } + + // Data representation + + Eigen::VectorXd serialize() const + { + /*! \f{equation*}{ () \rightarrow \mathbb{C}^{n \times 1} \f} + * + * Returns a serialized representation. + */ + + const ptrdiff_t dim = this->get_dimension(); + const Eigen::MatrixXcd A = this->get_ados_representation(); + + Eigen::VectorXd out = Eigen::VectorXd::Zero(this->get_dimension()); + ptrdiff_t kk = 0; + for (ptrdiff_t ii = 0; ii < this->shape; ii++) + { + for (ptrdiff_t jj = 0; jj < this->shape; jj++) + { + out(kk) = std::real(A(ii,jj)); + out(kk+1) = std::imag(A(ii,jj)); + kk = kk + 2; + } + } + + return out; + } + + void unserialize(const Eigen::VectorXd &vector) + { + /*! \f{equation*}{ (\mathbb{R}^{n \times 1}) \rightarrow () \f} + * + * Sets the GL object from a serialized vector. + */ + + const ptrdiff_t dim = vector.size(); + this->shape = static_cast(std::sqrt(dim/2)); + this->_data = Eigen::MatrixXcd::Zero(this->shape, this->shape); + + ptrdiff_t kk = 0; + for (ptrdiff_t ii = 0; ii < this->shape; ii++) + { + for (ptrdiff_t jj = 0; jj < this->shape; jj++) + { + this->_data(ii,jj) = std::complex(vector(kk), vector(kk+1)); + kk = kk + 2; + } + } + } + + std::complex operator()(const size_t index1, const size_t index2) const + { + /*! \f{equation*}{ (\mathbb{Z}, \mathbb{Z}) \rightarrow \mathbb{C} \f} + * + * Gets a value in the square matrix representation. + */ + + return _data(index1, index2); + } + + std::complex & operator()(const size_t index1, const size_t index2) + { + /*! \f{equation*}{ (\mathbb{Z}, \mathbb{Z}) \rightarrow \mathbb{C} \f} + * + * Gets a value in the square matrix representation. + */ + + return _data(index1, index2); + } + + size_t rows() const + { + // TODO: Delete this method + return shape; + } + + size_t cols() const + { + // TODO: Delete this method + return shape; + } + + GLC operator*(const GLC & other) const + { + /*! \f{equation*}{ (GLC, GLC) \rightarrow GLC \f} + * + * Group product. + */ + + assert(this->shape == other.shape); + Base out = _data * other._data; + return out; + } + + GLC & operator*=(const GLC & other) + { + /*! \f{equation*}{ (GLC, GLC) \rightarrow GLC \f} + * + * In place group product. + */ + + assert(this->shape == other.shape); + this->_data *= other._data; + return *this; + } + + friend std::ostream & operator<<(std::ostream& os, const GLC & other); + }; + + std::ostream & operator<<(std::ostream& os, const GLC & other) + { + /*! + * Overloads the "<<" stream insertion operator. + */ + + os << static_cast(other._data); + return os; + } + } +} + +#endif diff --git a/Lielab/domain/lieiii.hpp b/Lielab/domain/lieiii.hpp index 40a5ae0..3e0179a 100644 --- a/Lielab/domain/lieiii.hpp +++ b/Lielab/domain/lieiii.hpp @@ -13,21 +13,23 @@ namespace domain using std::conditional; using std::is_same; template -using lieiii = typename conditional::value, CN, - typename conditional::value, cn, - typename conditional::value, GL, - typename conditional::value, gl, - typename conditional::value, RN, - typename conditional::value, rn, - typename conditional::value, SE, - typename conditional::value, se, - typename conditional::value, SO, - typename conditional::value, so, - typename conditional::value, SP, - typename conditional::value, sp, - typename conditional::value, SU, - typename conditional::value, su, - Eigen::MatrixXcd>::type>::type>::type>::type>::type>::type>::type>::type>::type>::type>::type>::type>::type>::type; // Default false is complex matrix, throw compile time error instead?? +using lieiii = typename conditional::value, CN, + typename conditional::value, cn, + typename conditional::value, GL, + typename conditional::value, gl, + typename conditional::value, GLC, + typename conditional::value, glc, + typename conditional::value, RN, + typename conditional::value, rn, + typename conditional::value, SE, + typename conditional::value, se, + typename conditional::value, SO, + typename conditional::value, so, + typename conditional::value, SP, + typename conditional::value, sp, + typename conditional::value, SU, + typename conditional::value, su, + Eigen::MatrixXcd>::type>::type>::type>::type>::type>::type>::type>::type>::type>::type>::type>::type>::type>::type>::type>::type; // Default false is complex matrix, throw compile time error instead?? } } diff --git a/Lielab/functions/cayley1.hpp b/Lielab/functions/cayley1.hpp index cb6d324..1e5ac41 100644 --- a/Lielab/functions/cayley1.hpp +++ b/Lielab/functions/cayley1.hpp @@ -45,6 +45,21 @@ namespace Lielab return (Id + m/2.0)*(Id - m/2.0).inverse(); } + template <> + Lielab::domain::GLC cayley1(const Lielab::domain::glc & a) + { + /* + * Cayley1 overload for glc + * + * Needed since glc is complex + */ + + const Eigen::MatrixXcd m = a.get_ados_representation(); + const Eigen::MatrixXcd Id = Eigen::MatrixXcd::Identity(a.shape, a.shape); + + return (Id + m/2.0)*(Id - m/2.0).inverse(); + } + template <> Lielab::domain::SU cayley1(const Lielab::domain::su & a) { diff --git a/Lielab/topos/functions/Ad.hpp b/Lielab/topos/functions/Ad.hpp index 3afc2cc..7f1fc68 100644 --- a/Lielab/topos/functions/Ad.hpp +++ b/Lielab/topos/functions/Ad.hpp @@ -29,6 +29,11 @@ Lielab::domain::CompositeAlgebra Ad(const Lielab::domain::CompositeManifold & A, out.space.push_back(Lielab::functions::Ad(std::get(A.space[ii]), std::get(b.space[ii]))); } + else if (ind == Lielab::domain::CompositeAlgebra::INDEX_glc) + { + out.space.push_back(Lielab::functions::Ad(std::get(A.space[ii]), + std::get(b.space[ii]))); + } else if (ind == Lielab::domain::CompositeAlgebra::INDEX_rn) { out.space.push_back(Lielab::functions::Ad(std::get(A.space[ii]), diff --git a/Lielab/topos/functions/cayley1.hpp b/Lielab/topos/functions/cayley1.hpp index 9ef0648..31751c8 100644 --- a/Lielab/topos/functions/cayley1.hpp +++ b/Lielab/topos/functions/cayley1.hpp @@ -27,6 +27,10 @@ namespace Lielab { out.space.push_back(Lielab::functions::cayley1(std::get(la.space[ii]))); } + else if (ind == Lielab::domain::CompositeAlgebra::INDEX_glc) + { + out.space.push_back(Lielab::functions::cayley1(std::get(la.space[ii]))); + } else if (ind == Lielab::domain::CompositeAlgebra::INDEX_rn) { out.space.push_back(Lielab::functions::cayley1(std::get(la.space[ii]))); diff --git a/Lielab/topos/functions/dexp.hpp b/Lielab/topos/functions/dexp.hpp index 8727198..123b7b0 100644 --- a/Lielab/topos/functions/dexp.hpp +++ b/Lielab/topos/functions/dexp.hpp @@ -29,6 +29,11 @@ Lielab::domain::CompositeAlgebra dexp(const Lielab::domain::CompositeAlgebra & a out.space.push_back(Lielab::functions::dexp(std::get(a.space[ii]), std::get(b.space[ii]), order)); } + else if (ind == Lielab::domain::CompositeAlgebra::INDEX_glc) + { + out.space.push_back(Lielab::functions::dexp(std::get(a.space[ii]), + std::get(b.space[ii]), order)); + } else if (ind == Lielab::domain::CompositeAlgebra::INDEX_rn) { out.space.push_back(Lielab::functions::dexp(std::get(a.space[ii]), diff --git a/Lielab/topos/functions/exp.hpp b/Lielab/topos/functions/exp.hpp index ef45357..ca401b5 100644 --- a/Lielab/topos/functions/exp.hpp +++ b/Lielab/topos/functions/exp.hpp @@ -28,6 +28,10 @@ Lielab::domain::CompositeManifold exp(const Lielab::domain::CompositeAlgebra & l { out.space.push_back(Lielab::functions::exp(std::get(la.space[ii]))); } + else if (ind == Lielab::domain::CompositeAlgebra::INDEX_glc) + { + out.space.push_back(Lielab::functions::exp(std::get(la.space[ii]))); + } else if (ind == Lielab::domain::CompositeAlgebra::INDEX_rn) { out.space.push_back(Lielab::functions::exp(std::get(la.space[ii]))); diff --git a/Lielab/topos/functions/log.hpp b/Lielab/topos/functions/log.hpp index 255737f..e1985b5 100644 --- a/Lielab/topos/functions/log.hpp +++ b/Lielab/topos/functions/log.hpp @@ -28,6 +28,10 @@ Lielab::domain::CompositeAlgebra log(const Lielab::domain::CompositeManifold &M) { out.space.push_back(Lielab::functions::log(std::get(M.space[ii]))); } + else if (ind == Lielab::domain::CompositeManifold::INDEX_GLC) + { + out.space.push_back(Lielab::functions::log(std::get(M.space[ii]))); + } else if (ind == Lielab::domain::CompositeManifold::INDEX_RN) { out.space.push_back(Lielab::functions::log(std::get(M.space[ii]))); diff --git a/python/cppLielab.cpp b/python/cppLielab.cpp index e6317c3..8ff9aab 100644 --- a/python/cppLielab.cpp +++ b/python/cppLielab.cpp @@ -169,6 +169,67 @@ PYBIND11_MODULE(cppLielab, m) { { return matstr(self._data); }); + + + /*! + * Bindings for Lielab::domain::glc + */ + + py::class_(m_domain, "glc") + .def_readwrite("_data", &Lielab::domain::glc::_data) + .def_readonly_static("abelian", &Lielab::domain::glc::abelian) + .def(py::init<>()) + .def(py::init()) + .def(py::init()) + .def("basis", &Lielab::domain::glc::basis) + .def("project", &Lielab::domain::glc::project) + .def("get_dimension", &Lielab::domain::glc::get_dimension) + .def("get_vector", &Lielab::domain::glc::get_vector) + .def("get_ados_representation", &Lielab::domain::glc::get_ados_representation) + .def("set_vector", &Lielab::domain::glc::set_vector) + .def("__call__", [](Lielab::domain::glc & self, const size_t index) + { + return self(index); + }) + .def("__call__", [](const Lielab::domain::glc & self, const size_t index1, const size_t index2) + { + return self(index1, index2); + }) + .def(py::self + py::self) + .def(py::self += py::self) + .def(py::self - py::self) + .def(py::self -= py::self) + .def(-py::self) + .def(py::self * int()) + .def(py::self * double()) + // .def(py::self * std::complex()) // TODO: Complex integers seem bugged in Eigen right now + .def(py::self * std::complex()) + .def(py::self * py::self) + .def(py::self *= int()) + .def(py::self *= double()) + // .def(py::self *= std::complex()) // TODO: Complex integers seem bugged in Eigen right now + .def(py::self *= std::complex()) + .def(py::self *= py::self) + .def(py::self / int()) + .def(py::self / double()) + // .def(py::self / std::complex()) + .def(py::self / std::complex()) + .def(py::self /= int()) + .def(py::self /= double()) + // .def(py::self /= std::complex()) + .def(py::self /= std::complex()) + .def(int() * py::self) + .def(double() * py::self) + // .def(std::complex() * py::self) + .def(std::complex() * py::self) + .def("__repr__", [](const Lielab::domain::glc & self) + { + return ""; + }) + .def("__str__", [](const Lielab::domain::glc & self) + { + return matstr(self._data); + }); /*! * Bindings for Lielab::domain::rn @@ -496,6 +557,33 @@ PYBIND11_MODULE(cppLielab, m) { { return matstr(self._data); }); + + py::class_(m_domain, "GLC") + .def_readwrite("_data", &Lielab::domain::GLC::_data) + .def_readonly_static("abelian", &Lielab::domain::GLC::abelian) + .def(py::init<>()) + .def(py::init()) + .def(py::init()) + .def("project", &Lielab::domain::GLC::project) + .def("get_dimension", &Lielab::domain::GLC::get_dimension) + .def("get_ados_representation", &Lielab::domain::GLC::get_ados_representation) + .def("inverse", &Lielab::domain::GLC::inverse) + .def("serialize", &Lielab::domain::GLC::serialize) + .def("unserialize", &Lielab::domain::GLC::unserialize) + .def("__call__", [](const Lielab::domain::GLC & self, const size_t index1, const size_t index2) + { + return self(index1, index2); + }) + .def(py::self * py::self) + .def(py::self *= py::self) + .def("__repr__", [](const Lielab::domain::GLC & self) + { + return ""; + }) + .def("__str__", [](const Lielab::domain::GLC & self) + { + return matstr(self._data); + }); py::class_(m_domain, "RN") .def_readwrite("_data", &Lielab::domain::RN::_data) @@ -710,6 +798,7 @@ PYBIND11_MODULE(cppLielab, m) { py::module m_functions = m.def_submodule("functions", "The functions submodule."); m_functions.def("pair", &Lielab::functions::pair, "The pair function."); m_functions.def("pair", &Lielab::functions::pair, "The pair function."); + m_functions.def("pair", &Lielab::functions::pair, "The pair function."); m_functions.def("pair", &Lielab::functions::pair, "The pair function."); m_functions.def("pair", &Lielab::functions::pair, "The pair function."); m_functions.def("pair", &Lielab::functions::pair, "The pair function."); @@ -723,6 +812,7 @@ PYBIND11_MODULE(cppLielab, m) { m_functions.def("factorial", &Lielab::functions::factorial, "The factorial function"); m_functions.def("Ad", &Lielab::functions::Ad, "The Ad function."); m_functions.def("Ad", &Lielab::functions::Ad, "The Ad function."); + m_functions.def("Ad", &Lielab::functions::Ad, "The Ad function."); m_functions.def("Ad", &Lielab::functions::Ad, "The Ad function."); m_functions.def("Ad", &Lielab::functions::Ad, "The Ad function."); m_functions.def("Ad", &Lielab::functions::Ad, "The Ad function."); @@ -730,6 +820,7 @@ PYBIND11_MODULE(cppLielab, m) { m_functions.def("Ad", &Lielab::functions::Ad, "The Ad function."); m_functions.def("commutator", &Lielab::functions::commutator, "The commutator function."); m_functions.def("commutator", &Lielab::functions::commutator, "The commutator function."); + m_functions.def("commutator", &Lielab::functions::commutator, "The commutator function."); m_functions.def("commutator", &Lielab::functions::commutator, "The commutator function."); m_functions.def("commutator", &Lielab::functions::commutator, "The commutator function."); m_functions.def("commutator", &Lielab::functions::commutator, "The commutator function."); @@ -737,6 +828,7 @@ PYBIND11_MODULE(cppLielab, m) { m_functions.def("commutator", &Lielab::functions::commutator, "The commutator function."); // m_functions.def("cayley1", &Lielab::functions::cayley1, "The cayley1 function."); // m_functions.def("cayley1", &Lielab::functions::cayley1, "The cayley1 function."); + // m_functions.def("cayley1", &Lielab::functions::cayley1, "The cayley1 function."); m_functions.def("cayley1", &Lielab::functions::cayley1, "The cayley1 function."); // m_functions.def("cayley1", &Lielab::functions::cayley1, "The cayley1 function."); m_functions.def("cayley1", &Lielab::functions::cayley1, "The cayley1 function."); @@ -744,6 +836,7 @@ PYBIND11_MODULE(cppLielab, m) { m_functions.def("cayley1", &Lielab::functions::cayley1, "The cayley1 function."); // m_functions.def("cayley2", &Lielab::functions::cayley2, "The cayley2 function."); // m_functions.def("cayley2", &Lielab::functions::cayley2, "The cayley2 function."); + // m_functions.def("cayley2", &Lielab::functions::cayley2, "The cayley2 function."); m_functions.def("cayley2", &Lielab::functions::cayley2, "The cayley2 function."); // m_functions.def("cayley2", &Lielab::functions::cayley1, "The cayley2 function."); m_functions.def("cayley2", &Lielab::functions::cayley2, "The cayley2 function."); @@ -751,6 +844,7 @@ PYBIND11_MODULE(cppLielab, m) { // m_functions.def("cayley2", &Lielab::functions::cayley2, "The cayley2 function."); // m_functions.def("Killing", &Lielab::functions::Killing, "The Killing function."); // m_functions.def("Killing", &Lielab::functions::Killing, "The Killing function."); + // m_functions.def("Killing", &Lielab::functions::Killing, "The Killing function."); m_functions.def("Killing", &Lielab::functions::Killing, "The Killing function."); // m_functions.def("Killing", &Lielab::functions::Killing, "The Killing function."); m_functions.def("Killing", &Lielab::functions::Killing, "The Killing function."); @@ -758,6 +852,7 @@ PYBIND11_MODULE(cppLielab, m) { m_functions.def("Killing", &Lielab::functions::Killing, "The Killing function."); // TODO: Might be wrong m_functions.def("Killingform", &Lielab::functions::Killingform, "The Killingform function."); m_functions.def("Killingform", &Lielab::functions::Killingform, "The Killingform function."); + m_functions.def("Killingform", &Lielab::functions::Killingform, "The Killingform function."); // m_functions.def("Killingform", &Lielab::functions::Killingform, "The Killingform function."); m_functions.def("Killingform", &Lielab::functions::Killingform, "The Killingform function."); m_functions.def("Killingform", &Lielab::functions::Killingform, "The Killingform function."); @@ -765,6 +860,7 @@ PYBIND11_MODULE(cppLielab, m) { m_functions.def("Killingform", &Lielab::functions::Killingform, "The Killingform function."); // TODO: Might be wrong // m_functions.def("ad", &Lielab::functions::ad, "The ad function."); // m_functions.def("ad", &Lielab::functions::ad, "The ad function."); + // m_functions.def("ad", &Lielab::functions::ad, "The ad function."); m_functions.def("ad", &Lielab::functions::ad, "The ad function."); // m_functions.def("ad", &Lielab::functions::ad, "The ad function."); m_functions.def("ad", &Lielab::functions::ad, "The ad function."); @@ -782,6 +878,7 @@ PYBIND11_MODULE(cppLielab, m) { // m_functions.def("coad", &Lielab::functions::coad, "The coad function."); // TODO: su needs basis() m_functions.def("exp", &Lielab::functions::exp, "The exponential function."); m_functions.def("exp", &Lielab::functions::exp, "The exponential function."); + m_functions.def("exp", &Lielab::functions::exp, "The exponential function."); m_functions.def("exp", &Lielab::functions::exp, "The exponential function."); m_functions.def("exp", &Lielab::functions::exp, "The exponential function."); m_functions.def("exp", &Lielab::functions::exp, "The exponential function."); @@ -789,6 +886,7 @@ PYBIND11_MODULE(cppLielab, m) { m_functions.def("exp", &Lielab::functions::exp, "The exponential function."); m_functions.def("log", &Lielab::functions::log, "The logarithm function.", py::arg("G"), py::arg("optimize") = false); m_functions.def("log", &Lielab::functions::log, "The logarithm function.", py::arg("G"), py::arg("optimize") = false); + m_functions.def("log", &Lielab::functions::log, "The logarithm function.", py::arg("G"), py::arg("optimize") = false); m_functions.def("log", &Lielab::functions::log, "The logarithm function.", py::arg("G"), py::arg("optimize") = false); m_functions.def("log", &Lielab::functions::log, "The logarithm function.", py::arg("G"), py::arg("optimize") = false); m_functions.def("log", &Lielab::functions::log, "The logarithm function.", py::arg("G"), py::arg("optimize") = false); @@ -797,6 +895,7 @@ PYBIND11_MODULE(cppLielab, m) { m_functions.def("bernoulli", &Lielab::functions::bernoulli, "The bernoulli function."); // m_functions.def("dcayley1inv", &Lielab::functions::dcayley1inv, "The dcayley1inv function."); m_functions.def("dcayley1inv", &Lielab::functions::dcayley1inv, "The dcayley1inv function."); + // m_functions.def("dcayley1inv", &Lielab::functions::dcayley1inv, "The dcayley1inv function."); m_functions.def("dcayley1inv", &Lielab::functions::dcayley1inv, "The dcayley1inv function."); // m_functions.def("dcayley1inv", &Lielab::functions::dcayley1inv, "The dcayley1inv function."); m_functions.def("dcayley1inv", &Lielab::functions::dcayley1inv, "The dcayley1inv function."); @@ -804,6 +903,7 @@ PYBIND11_MODULE(cppLielab, m) { m_functions.def("dcayley1inv", &Lielab::functions::dcayley1inv, "The dcayley1inv function."); m_functions.def("dexp", &Lielab::functions::dexp, "The dexp function.", py::arg("a"), py::arg("b"), py::arg("order") = 5); m_functions.def("dexp", &Lielab::functions::dexp, "The dexp function.", py::arg("a"), py::arg("b"), py::arg("order") = 5); + m_functions.def("dexp", &Lielab::functions::dexp, "The dexp function.", py::arg("a"), py::arg("b"), py::arg("order") = 5); m_functions.def("dexp", &Lielab::functions::dexp, "The dexp function.", py::arg("a"), py::arg("b"), py::arg("order") = 5); m_functions.def("dexp", &Lielab::functions::dexp, "The dexp function.", py::arg("a"), py::arg("b"), py::arg("order") = 5); m_functions.def("dexp", &Lielab::functions::dexp, "The dexp function.", py::arg("a"), py::arg("b"), py::arg("order") = 5); @@ -811,6 +911,7 @@ PYBIND11_MODULE(cppLielab, m) { m_functions.def("dexp", &Lielab::functions::dexp, "The dexp function.", py::arg("a"), py::arg("b"), py::arg("order") = 5); m_functions.def("dexpinv", &Lielab::functions::dexpinv, "The dexpinv function.", py::arg("a"), py::arg("b"), py::arg("order") = 5); m_functions.def("dexpinv", &Lielab::functions::dexpinv, "The dexpinv function.", py::arg("a"), py::arg("b"), py::arg("order") = 5); + m_functions.def("dexpinv", &Lielab::functions::dexpinv, "The dexpinv function.", py::arg("a"), py::arg("b"), py::arg("order") = 5); m_functions.def("dexpinv", &Lielab::functions::dexpinv, "The dexpinv function.", py::arg("a"), py::arg("b"), py::arg("order") = 5); m_functions.def("dexpinv", &Lielab::functions::dexpinv, "The dexpinv function.", py::arg("a"), py::arg("b"), py::arg("order") = 5); m_functions.def("dexpinv", &Lielab::functions::dexpinv, "The dexpinv function.", py::arg("a"), py::arg("b"), py::arg("order") = 5); diff --git a/python/tests/domain/tests/test_domain.py b/python/tests/domain/tests/test_domain.py index 87f7170..d257929 100644 --- a/python/tests/domain/tests/test_domain.py +++ b/python/tests/domain/tests/test_domain.py @@ -127,6 +127,66 @@ def test_gl(): assert_domain(commutator(x,y), -commutator(y,x)) +def test_glc(): + """ + Checks the basic glc class + """ + + from lielab.domain import glc + from lielab.functions import commutator + + one = glc(1) + two = glc(2) + three = glc(3) + four = glc(4) + five = glc(5) + six = glc(6) + seven = glc(7) + eight = glc(8) + + # Dimensions + assert one.get_dimension() == 2 + assert two.get_dimension() == 8 + assert three.get_dimension() == 18 + assert four.get_dimension() == 32 + assert five.get_dimension() == 50 + assert six.get_dimension() == 72 + assert seven.get_dimension() == 98 + assert eight.get_dimension() == 128 + + x = glc.basis(0,4) + y = glc.basis(1,4) + z = glc.basis(2,4) + zero = x*0 + + # Unary subtraction + -x + + # Vector Addition + x + y + assert_domain(x + y, y + x) + + # Vector Subtraction + x - y + assert_domain(x - y, -(y - x)) + + # Scalar Multiplication + an_int*x + a_double*x + x*an_int + x*a_double + assert_domain(an_int*x, x*an_int) + assert_domain(a_double*x, x*a_double) + + # Scalar division + x/an_int + x/a_double + + # Vector multiplication + x * y + assert_domain(commutator(x,y), -commutator(y,x)) + + def test_rn(): """ Checks the basic rn class @@ -805,6 +865,26 @@ def test_GL(): x.inverse() +def test_GLC(): + """ + Tests GLC against well-known identities. + """ + + from lielab.domain import glc + from lielab.functions import exp + + x = exp(glc.basis(0,2)) + y = exp(glc.basis(1,2)) + z = exp(glc.basis(2,2)) + + # Group multiplication + x*y + assert_domain(x*y, (y.inverse()*x.inverse()).inverse()) + + # Group inverse + x.inverse() + + def test_RN(): """ Tests RN against well-known identities. diff --git a/tests/test_domain.cpp b/tests/test_domain.cpp index d7d0e0d..9b2ebb1 100644 --- a/tests/test_domain.cpp +++ b/tests/test_domain.cpp @@ -198,6 +198,42 @@ TEST_CASE("gl algebra", "[domain]") CHECK(eight.get_dimension() == 64); } + +TEST_CASE("glc algebra", "[domain]") +{ + /*! + * Tests the glc algebra. + */ + + for (size_t shape = 1; shape <= TEST_UP_TO_THIS_SHAPE; shape++) + { + const size_t D = Lielab::domain::glc::basis(0, shape).get_dimension(); + + // Construct the glc basis + std::vector basis; + for (size_t ii = 0; ii < D; ii++) + { + basis.push_back(Lielab::domain::glc::basis(ii, shape)); + } + + is_algebra(basis); + is_liealgebra(basis); + } + + Lielab::domain::glc one(1), two(2), three(3), four(4), five(5), six(6), seven(7), eight(8); + + // Dimensions + CHECK(one.get_dimension() == 2); + CHECK(two.get_dimension() == 8); + CHECK(three.get_dimension() == 18); + CHECK(four.get_dimension() == 32); + CHECK(five.get_dimension() == 50); + CHECK(six.get_dimension() == 72); + CHECK(seven.get_dimension() == 98); + CHECK(eight.get_dimension() == 128); +} + + TEST_CASE("rn algebra", "[domain]") { /*! @@ -519,6 +555,29 @@ TEST_CASE("GL", "[domain]") } } +TEST_CASE("GLC", "[domain]") +{ + /*! + * Tests GLC against well-known identities. + */ + + for (size_t shape = 1; shape <= TEST_UP_TO_THIS_SHAPE; shape++) + { + const size_t D = Lielab::domain::glc::basis(0, shape).get_dimension(); + + // Construct the GL elements + std::vector elements; + for (size_t ii = 0; ii < D; ii++) + { + elements.push_back(Lielab::functions::exp(Lielab::domain::glc::basis(ii, shape))); + } + + const Lielab::domain::GLC identity(shape); + + is_group(elements, identity); + } +} + TEST_CASE("RN", "[domain]") { /*!