1 #ifndef LOCALREGRESSION_H_ 2 #define LOCALREGRESSION_H_ 3 4 #include "MUQ/config.h" 5 6 #if MUQ_HAS_PARCER 7 #include <parcer/Communicator.h> 8 #endif 9 10 #include "MUQ/Modeling/ModPiece.h" 11 #include "MUQ/Modeling/Flann/FlannCache.h" 12 13 #include "MUQ/Approximation/Regression/Regression.h" 14 15 namespace muq { 16 namespace Approximation { 17 class LocalRegression : public muq::Modeling::ModPiece { 18 public: 19 20 /** 21 @param[in] function The function we wish to approximate with a local polynomial 22 @param[in] pt Options for the regression 23 */ 24 LocalRegression(std::shared_ptr<muq::Modeling::ModPiece> function, boost::property_tree::ptree& pt); 25 26 #if MUQ_HAS_PARCER 27 /** 28 @param[in] function The function we wish to approximate with a local polynomial 29 @param[in] pt Options for the regression 30 @param[in] comm The parcer communicator 31 */ 32 LocalRegression(std::shared_ptr<muq::Modeling::ModPiece> function, boost::property_tree::ptree& pt, std::shared_ptr<parcer::Communicator> comm); 33 #endif 34 35 #if MUQ_HAS_PARCER 36 ~LocalRegression(); 37 #else 38 ~LocalRegression() = default; 39 #endif 40 41 /// Add some points to the cache 42 /** 43 @param[in] inputs Points to add to the cache 44 */ 45 void Add(std::vector<Eigen::VectorXd> const& inputs) const; 46 47 /// Add a single point to the cache 48 /** 49 @param[in] input Point to add to the cache 50 \return The function result at the new point 51 */ 52 Eigen::VectorXd Add(Eigen::VectorXd const& input) const; 53 54 /// Get the total size of the cache 55 /** 56 \return The cache size 57 */ 58 unsigned int CacheSize() const; 59 60 /// Is a point in the cache? 61 /** 62 @param[in] point We want to know if this point is in the cache 63 \return true: it is in the cache, false: it is not in the cache 64 */ 65 bool InCache(Eigen::VectorXd const& point) const; 66 67 /// A point in the cache 68 /** 69 @param[in] index The index of the input 70 \return The input point in the cache associated with the index 71 */ 72 Eigen::VectorXd CachePoint(unsigned int const index) const; 73 74 /// Get the poisedness constant 75 /** 76 Get the poisedness constant associated with the nearest neighbors of the input point. 77 @param[in] input The input point 78 \return first: the point where the Lagrange polynomials are maximized, second: the poisedness constant, third: the index of the nearest neighbor closest to the new point 79 */ 80 std::tuple<Eigen::VectorXd, double, unsigned int> PoisednessConstant(Eigen::VectorXd const& input) const; 81 82 /// Get the poisedness constant 83 /** 84 Get the poisedness constant associated with the nearest neighbors of the input point. 85 @param[in] input The input point 86 @param[in] neighbors The nearest neighbors 87 \return first: the point where the Lagrange polynomials are maximized, second: the poisedness constant, third: the index of the nearest neighbor closest to the new point 88 */ 89 std::tuple<Eigen::VectorXd, double, unsigned int> PoisednessConstant(Eigen::VectorXd const& input, std::vector<Eigen::VectorXd> const& neighbors) const; 90 91 /// Get the error indicator 92 /** 93 Get the error indicator \f$\Lambda \sqrt{k} \Delta^{p+1}\f$ 94 @param[in] input The input point 95 \return first: the error indicator, second: the radius of the ball 96 */ 97 std::pair<double, double> ErrorIndicator(Eigen::VectorXd const& input) const; 98 99 /// Get the error indicator 100 /** 101 Get the error indicator \f$\Lambda \sqrt{k} \Delta^{p+1}\f$ 102 @param[in] input The input point 103 @param[in] neighbors The nearest neighbors 104 \return first: the error indicator, second: the radius of the ball 105 */ 106 std::pair<double, double> ErrorIndicator(Eigen::VectorXd const& input, std::vector<Eigen::VectorXd> const& neighbors) const; 107 108 /// Get the number of nearest neighbors 109 /** 110 @param[in] input We want the \f$k\f$ nearest neighbors to this point 111 @param[out] neighbors The \f$k\f$ nearest neighbors 112 */ 113 void NearestNeighbors(Eigen::VectorXd const& input, std::vector<Eigen::VectorXd>& neighbors) const; 114 115 /// Get the number of nearest neighbors (with result) 116 /** 117 @param[in] input We want the \f$k\f$ nearest neighbors to this point 118 @param[out] neighbors The \f$k\f$ nearest neighbors 119 @param[out] results The corresponding output of the function 120 */ 121 void NearestNeighbors(Eigen::VectorXd const& input, std::vector<Eigen::VectorXd>& neighbors, std::vector<Eigen::VectorXd>& result) const; 122 123 Eigen::VectorXd EvaluateRegressor(Eigen::VectorXd const& input, std::vector<Eigen::VectorXd> const& neighbors, std::vector<Eigen::VectorXd> const& result) const; 124 125 #if MUQ_HAS_PARCER 126 void Probe() const; 127 #endif 128 129 /// The number of nearest neighbors to use by the regressor 130 const unsigned int kn; 131 132 /// Get the centroid of the cache 133 /** 134 \return The cache centroid 135 */ 136 Eigen::VectorXd CacheCentroid() const; 137 138 /// Clear the cache 139 void ClearCache(); 140 141 /// Polynomial order 142 unsigned int Order() const; 143 144 private: 145 146 virtual void EvaluateImpl(muq::Modeling::ref_vector<Eigen::VectorXd> const& inputs) override; 147 148 /// Fit the regression to the nearest neighbors 149 void FitRegression(Eigen::VectorXd const& input) const; 150 151 /// Set up the regressor 152 void SetUp(std::shared_ptr<muq::Modeling::ModPiece> function, boost::property_tree::ptree& pt); 153 154 /// A cache containing previous model evaluations 155 std::shared_ptr<muq::Modeling::FlannCache> cache; 156 157 /// A regressor 158 std::shared_ptr<Regression> reg; 159 160 #if MUQ_HAS_PARCER 161 std::shared_ptr<parcer::Communicator> comm = nullptr; 162 const int tagSingle = 0; 163 #endif 164 }; 165 } // namespace Approximation 166 } // namespace muq 167 168 #endif 169