HybZono.hpp

Go to the documentation of this file.

#ifndef ZONOOPT_HYBZONO_HPP_
#define ZONOOPT_HYBZONO_HPP_
 
#include "SparseMatrixUtilities.hpp"
#include "BranchAndBound.hpp"
#include "BnbDataStructures.hpp"
#include "IntervalMatrix.hpp"
 
#include <cassert>
#include <limits>
#include <set>
 
namespace ZonoOpt
{
 
using namespace detail;
 
// forward declarations
class Zono;
class ConZono;
 
class HybZono
{
    public:
        
        // constructors
 
        HybZono() = default;
 
        HybZono(const Eigen::SparseMatrix<zono_float>& Gc, const Eigen::SparseMatrix<zono_float>& Gb, const Eigen::Vector<zono_float, -1>& c,
            const Eigen::SparseMatrix<zono_float>& Ac, const Eigen::SparseMatrix<zono_float>& Ab, const Eigen::Vector<zono_float, -1>& b,
            const bool zero_one_form=false, const bool sharp=false);
 
        // virtual destructor
        virtual ~HybZono() = default;
 
        void set(const Eigen::SparseMatrix<zono_float>& Gc, const Eigen::SparseMatrix<zono_float>& Gb, const Eigen::Vector<zono_float, -1>& c,
            const Eigen::SparseMatrix<zono_float>& Ac, const Eigen::SparseMatrix<zono_float>& Ab, const Eigen::Vector<zono_float, -1>& b,
            bool zero_one_form=false, bool sharp=false);
 
        virtual HybZono* clone() const;
 
        // get methods
 
        virtual int get_n() const { return this->n; }
 
        virtual int get_nC() const { return this->nC; }
 
        virtual int get_nG() const { return this->nG; }
 
        virtual int get_nGc() const { return this->nGc; }
 
        virtual int get_nGb() const { return this->nGb; }
 
        virtual Eigen::SparseMatrix<zono_float> get_Gc() const { return this->Gc; }
 
        virtual Eigen::SparseMatrix<zono_float> get_Gb() const { return this->Gb; }
 
        virtual Eigen::SparseMatrix<zono_float> get_G() const { return this->G; }
 
        virtual Eigen::SparseMatrix<zono_float> get_Ac() const { return this->Ac; }
 
        virtual Eigen::SparseMatrix<zono_float> get_Ab() const { return this->Ab; }
 
        virtual Eigen::SparseMatrix<zono_float> get_A() const { return this->A; }
 
        virtual Eigen::Vector<zono_float, -1> get_c() const { return this->c; }
        
        virtual Eigen::Vector<zono_float, -1> get_b() const { return this->b; }
 
        virtual bool is_0_1_form() const { return this->zero_one_form; }
 
        bool is_sharp() const { return this->sharp; }
 
        virtual void convert_form();
 
        virtual std::unique_ptr<HybZono> remove_redundancy(int contractor_iter=10) const;
 
        virtual std::unique_ptr<ConZono> convex_relaxation() const;
 
        virtual std::unique_ptr<HybZono> complement(const zono_float delta_m = 100, const bool remove_redundancy=true, const OptSettings &settings=OptSettings(),
            std::shared_ptr<OptSolution>* solution=nullptr, const int n_leaves = std::numeric_limits<int>::max(), const int contractor_iter=10)
        {
            return do_complement(delta_m, remove_redundancy, settings, solution, n_leaves, contractor_iter);
        }
 
        // type checking
        bool is_point() const;
 
        bool is_zono() const;
 
        bool is_conzono() const;
 
        bool is_hybzono() const;
 
        bool is_empty_set() const;
        
        // display methods
        virtual std::string print() const;
 
        friend std::ostream& operator<<(std::ostream& os, const HybZono& Z);
 
        // optimization
        Eigen::Vector<zono_float, -1> optimize_over(
            const Eigen::SparseMatrix<zono_float> &P, const Eigen::Vector<zono_float, -1> &q, zono_float c=0,
            const OptSettings &settings=OptSettings(), std::shared_ptr<OptSolution>* solution=nullptr,
            const WarmStartParams& warm_start_params=WarmStartParams()) const
        {
            return do_optimize_over(P, q, c, settings, solution, warm_start_params);
        }
 
        Eigen::Vector<zono_float, -1> project_point(const Eigen::Vector<zono_float, -1>& x,
            const OptSettings &settings=OptSettings(), std::shared_ptr<OptSolution>* solution=nullptr,
            const WarmStartParams& warm_start_params=WarmStartParams()) const
        {
            return do_project_point(x, settings, solution, warm_start_params);
        }
 
        bool is_empty(const OptSettings &settings=OptSettings(),
            std::shared_ptr<OptSolution>* solution=nullptr, const WarmStartParams& warm_start_params=WarmStartParams()) const
        {
            return do_is_empty(settings, solution, warm_start_params);
        }
 
        zono_float support(const Eigen::Vector<zono_float, -1>& d, const OptSettings &settings=OptSettings(),
            std::shared_ptr<OptSolution>* solution=nullptr, const WarmStartParams& warm_start_params=WarmStartParams())
        {
            return do_support(d, settings, solution, warm_start_params);
        }
 
        bool contains_point(const Eigen::Vector<zono_float, -1>& x, const OptSettings &settings=OptSettings(),
            std::shared_ptr<OptSolution>* solution=nullptr, const WarmStartParams& warm_start_params=WarmStartParams()) const
        {
            return do_contains_point(x, settings, solution, warm_start_params);
        }
 
        Box bounding_box(const OptSettings &settings=OptSettings(), std::shared_ptr<OptSolution>* solution=nullptr,
            const WarmStartParams& warm_start_params=WarmStartParams())
        {
            return do_bounding_box(settings, solution, warm_start_params);
        }
 
        std::vector<std::unique_ptr<ConZono>> get_leaves(bool remove_redundancy=false, const OptSettings &settings=OptSettings(),
            std::shared_ptr<OptSolution>* solution=nullptr, int n_leaves = std::numeric_limits<int>::max(), int contractor_iter=10) const;
 
        // friend function declarations
        friend std::unique_ptr<HybZono> affine_map(const HybZono& Z,
            const Eigen::SparseMatrix<zono_float>& R, const Eigen::Vector<zono_float, -1>& s);
        friend std::unique_ptr<HybZono> affine_inclusion(const HybZono& Z, const IntervalMatrix& R, const Eigen::Vector<zono_float, -1>& s);
        friend std::unique_ptr<HybZono> project_onto_dims(const HybZono& Z, const std::vector<int>& dims);
        friend std::unique_ptr<HybZono> minkowski_sum(const HybZono& Z1, HybZono& Z2);
        friend std::unique_ptr<HybZono> pontry_diff(HybZono& Z1, Zono& Z2, bool exact);
        friend std::unique_ptr<HybZono> intersection(const HybZono& Z1, HybZono& Z2, 
            const Eigen::SparseMatrix<zono_float>& R);
        friend std::unique_ptr<HybZono> intersection_over_dims(const HybZono& Z1, HybZono& Z2, 
            const std::vector<int>& dims);
        friend std::unique_ptr<HybZono> halfspace_intersection(HybZono& Z, const Eigen::SparseMatrix<zono_float>& H, 
            const Eigen::Vector<zono_float, -1>& f, const Eigen::SparseMatrix<zono_float>& R);
        friend std::unique_ptr<HybZono> union_of_many(const std::vector<std::shared_ptr<HybZono>>& Zs, bool preserve_sharpness, bool expose_indicators);
        friend std::unique_ptr<ConZono> convex_hull(const std::vector<std::shared_ptr<HybZono>>& Zs);
        friend std::unique_ptr<HybZono> cartesian_product(const HybZono& Z1, HybZono& Z2);
        friend std::unique_ptr<HybZono> constrain(HybZono& Z, const Eigen::SparseMatrix<zono_float>& H,
            const Eigen::Vector<zono_float, -1>& f, char direction, const Eigen::SparseMatrix<zono_float>& R);
        friend std::unique_ptr<HybZono> set_diff(const HybZono& Z1, HybZono& Z2, zono_float delta_m, bool remove_redundancy,
            const OptSettings &settings, std::shared_ptr<OptSolution>* solution, const WarmStartParams& warm_start_params, int n_leaves, int contractor_iter);
        friend std::unique_ptr<HybZono> vrep_2_hybzono(const std::vector<Eigen::Matrix<zono_float, -1, -1>> &Vpolys, bool expose_indicators);
        friend std::unique_ptr<HybZono> zono_union_2_hybzono(std::vector<std::shared_ptr<Zono>> &Zs, bool expose_indicators);
 
        // operator overloading
 
        std::unique_ptr<HybZono> operator+(HybZono& other) const;
 
        std::unique_ptr<HybZono> operator+(const Eigen::Vector<zono_float, -1>& v) const;
 
        friend std::unique_ptr<HybZono> operator+(const Eigen::Vector<zono_float, -1>& v, HybZono& Z);
 
        std::unique_ptr<HybZono> operator+(const Box& box) const;
 
        friend std::unique_ptr<HybZono> operator+(const Box& box, HybZono& Z);
 
        void operator+=(HybZono& other);
 
        void operator+=(const Eigen::Vector<zono_float, -1>& v);
 
        void operator+=(const Box& box);
 
        friend std::unique_ptr<HybZono> operator*(const Eigen::SparseMatrix<zono_float>& R, const HybZono& Z);
 
        friend std::unique_ptr<HybZono> operator*(const Eigen::Matrix<zono_float, -1, -1>& R, const HybZono& Z);
 
        friend std::unique_ptr<HybZono> operator*(const IntervalMatrix& R, const HybZono& Z);
 
        std::unique_ptr<HybZono> operator*(zono_float f) const;
 
        friend std::unique_ptr<HybZono> operator*(zono_float f, const HybZono& Z);
 
        void operator*=(zono_float f);
 
        std::unique_ptr<HybZono> operator-(Zono& other);
 
        std::unique_ptr<HybZono> operator-(const Eigen::Vector<zono_float, -1>& v);
 
        std::unique_ptr<HybZono> operator-(const Box& box);
 
        void operator-=(Zono& other);
 
        void operator-=(const Eigen::Vector<zono_float, -1>& v);
 
        void operator-=(const Box& box);
 
        std::unique_ptr<HybZono> operator*(HybZono& other) const;
 
        std::unique_ptr<HybZono> operator*(const Box& box) const;
 
        friend std::unique_ptr<HybZono> operator*(const Box& box, HybZono& Z);
 
        void operator*=(HybZono& other);
 
        void operator*=(const Box& box);
 
        std::unique_ptr<HybZono> operator&(HybZono& other) const;
 
        std::unique_ptr<HybZono> operator|(HybZono& other) const;
 
        std::unique_ptr<HybZono> operator-() const; // unary minus
 
    protected:
 
        // fields
 
        Eigen::SparseMatrix<zono_float> G = Eigen::SparseMatrix<zono_float>(0, 0);
 
        Eigen::SparseMatrix<zono_float> Gc = Eigen::SparseMatrix<zono_float>(0, 0);
 
        Eigen::SparseMatrix<zono_float> Gb = Eigen::SparseMatrix<zono_float>(0, 0);
 
        Eigen::SparseMatrix<zono_float> A = Eigen::SparseMatrix<zono_float>(0, 0);
 
        Eigen::SparseMatrix<zono_float> Ac = Eigen::SparseMatrix<zono_float>(0, 0);
 
        Eigen::SparseMatrix<zono_float> Ab = Eigen::SparseMatrix<zono_float>(0, 0);
 
        Eigen::Vector<zono_float, -1> c = Eigen::Vector<zono_float, -1>(0);
 
        Eigen::Vector<zono_float, -1> b = Eigen::Vector<zono_float, -1>(0);
 
        int n = 0;
 
        int nG = 0;
 
        int nGc = 0;
 
        int nGb = 0;
 
        int nC = 0;
 
        bool zero_one_form = false;
 
        bool sharp = false;
 
        // methods
        virtual Eigen::Vector<zono_float, -1> do_optimize_over(
            const Eigen::SparseMatrix<zono_float> &P, const Eigen::Vector<zono_float, -1> &q, zono_float c,
            const OptSettings &settings, std::shared_ptr<OptSolution>* solution,
            const WarmStartParams& warm_start_params) const;
 
        virtual Eigen::Vector<zono_float, -1> do_project_point(const Eigen::Vector<zono_float, -1>& x,
            const OptSettings &settings, std::shared_ptr<OptSolution>* solution, const WarmStartParams& warm_start_params) const;
 
        virtual bool do_is_empty(const OptSettings &settings, std::shared_ptr<OptSolution>* solution,
            const WarmStartParams& warm_start_params) const;
 
        virtual zono_float do_support(const Eigen::Vector<zono_float, -1>& d, const OptSettings &settings,
            std::shared_ptr<OptSolution>* solution, const WarmStartParams& warm_start_params);
 
        virtual bool do_contains_point(const Eigen::Vector<zono_float, -1>& x, const OptSettings &settings,
            std::shared_ptr<OptSolution>* solution, const WarmStartParams& warm_start_params) const;
 
        virtual Box do_bounding_box(const OptSettings &settings, std::shared_ptr<OptSolution>* solution,
            const WarmStartParams& warm_start_params);
 
        virtual std::unique_ptr<HybZono> do_complement(zono_float, bool remove_redundancy, const OptSettings &settings,
            std::shared_ptr<OptSolution>* solution, int n_leaves, int contractor_iter);
 
 
        static void remove_cols(Eigen::SparseMatrix<zono_float>& M, const std::set<int>& idx_to_remove);
        static std::set<int> find_unused_generators(const Eigen::SparseMatrix<zono_float>& G, const Eigen::SparseMatrix<zono_float>& A);
        OptSolution mi_opt(const Eigen::SparseMatrix<zono_float>& P, const Eigen::Vector<zono_float, -1>& q,
            zono_float c, const Eigen::SparseMatrix<zono_float>& A, const Eigen::Vector<zono_float, -1>& b,
            const OptSettings &settings=OptSettings(), std::shared_ptr<OptSolution>* solution=nullptr,
            const WarmStartParams& warm_start_params=WarmStartParams()) const;
        std::vector<OptSolution> mi_opt_multisol(const Eigen::SparseMatrix<zono_float>& P, const Eigen::Vector<zono_float, -1>& q,
            zono_float c, const Eigen::SparseMatrix<zono_float>& A, const Eigen::Vector<zono_float, -1>& b, int n_sols,
            const OptSettings &settings=OptSettings(), std::shared_ptr<OptSolution>* solution=nullptr) const;
 
    private:
 
        void make_G_A();
        void set_Ac_Ab_from_A();
        std::vector<Eigen::Vector<zono_float, -1>> get_bin_leaves(const OptSettings &settings=OptSettings(), std::shared_ptr<OptSolution>* solution=nullptr,
            int n_leaves = std::numeric_limits<int>::max()) const;
        std::vector<std::pair<int, int>> get_simplifiable_constraints() const;
        void apply_constraint_simplification(const std::vector<std::pair<int, int>>& cons, Box& box);
        bool rescale_generators(MI_Box& box); // returns false if empty set detected
        void remove_generators(const std::set<int>& idx_c, const std::set<int>& idx_b, MI_Box& box);
        void remove_fixed_vars(MI_Box& box);
 
};
 
// forward delcarations
std::unique_ptr<HybZono> affine_map(const HybZono& Z,
    const Eigen::SparseMatrix<zono_float>& R, const Eigen::Vector<zono_float, -1>& s = Eigen::Vector<zono_float, -1>());
 
std::unique_ptr<HybZono> affine_inclusion(const HybZono& Z, const IntervalMatrix& R,
    const Eigen::Vector<zono_float, -1>& s=Eigen::Vector<zono_float, -1>());
 
std::unique_ptr<HybZono> project_onto_dims(const HybZono& Z, const std::vector<int>& dims);
 
std::unique_ptr<HybZono> minkowski_sum(const HybZono& Z1, HybZono& Z2);
 
std::unique_ptr<HybZono> pontry_diff(HybZono& Z1, Zono& Z2, bool exact=true);
 
std::unique_ptr<HybZono> intersection(const HybZono& Z1, HybZono& Z2,
    const Eigen::SparseMatrix<zono_float>& R=Eigen::SparseMatrix<zono_float>());
 
std::unique_ptr<HybZono> intersection_over_dims(const HybZono& Z1, HybZono& Z2,
    const std::vector<int>& dims);
 
std::unique_ptr<HybZono> halfspace_intersection(HybZono& Z, const Eigen::SparseMatrix<zono_float>& H,
    const Eigen::Vector<zono_float, -1>& f, const Eigen::SparseMatrix<zono_float>& R=Eigen::SparseMatrix<zono_float>());
 
std::unique_ptr<HybZono> union_of_many(const std::vector<std::shared_ptr<HybZono>>& Zs, bool preserve_sharpness=false, bool expose_indicators=false);
 
std::unique_ptr<ConZono> convex_hull(const std::vector<std::shared_ptr<HybZono>>& Zs);
 
std::unique_ptr<HybZono> cartesian_product(const HybZono& Z1, HybZono& Z2);
 
std::unique_ptr<HybZono> constrain(HybZono& Z, const Eigen::SparseMatrix<zono_float>& H,
            const Eigen::Vector<zono_float, -1>& f, char direction, const Eigen::SparseMatrix<zono_float>& R=Eigen::SparseMatrix<zono_float>());
 
std::unique_ptr<HybZono> set_diff(const HybZono& Z1, HybZono& Z2, zono_float delta_m = 100, bool remove_redundancy=true,
    const OptSettings &settings=OptSettings(), std::shared_ptr<OptSolution>* solution=nullptr,
    int n_leaves = std::numeric_limits<int>::max(), int contractor_iter = 10);
 
 
std::unique_ptr<HybZono> vrep_2_hybzono(const std::vector<Eigen::Matrix<zono_float, -1, -1>> &Vpolys, bool expose_indicators=false);
 
 
std::unique_ptr<HybZono> zono_union_2_hybzono(std::vector<std::shared_ptr<Zono>> &Zs, bool expose_indicators=false);
 
} // end namespace ZonoOpt
 
#endif